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Abstract 

This thesis presents a novel design methodology for the synthesis of automatic 
controllers, together with a computational environment—the Control Engineer’s 
Workbench—integrating a suite of programs that automatically analyze and de¬ 
sign controllers for high-performance, global control of nonlinear systems. This 
work demonstrates that difficult control synthesis tasks can be automated, using 
programs that actively exploit and efficiently represent knowledge of nonlinear dy¬ 
namics and phase space and effectively use the representation to guide and perform 
the control design. The Control Engineer’s Workbench combines powerful numer¬ 
ical and symbolic computations with spatial reasoning techniques. The two major 
programs in the Workbench—Phase Space Navigator and MAPS—work together 
to model and reason about the phase-space geometry and topology of a given sys¬ 
tem, to plan global control reference trajectories, and to navigate the system along 
the planned trajectories. They use a novel technique of “flow pipes” to group in¬ 
finite numbers of distinct behaviors into a manageable discrete set that becomes 
the basis for establishing the reference trajectories. 

As a demonstration of this approach, I exhibit the automatic design of a non¬ 
linear controller for a magnetic levitation system. The control system synthesized 
by the Workbench can stabilize a maglev vehicle with large initial displacements 
from an equilibrium and outperform the classical linear feedback design for the 
same system by a factor of 20. 
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Chapter 1 
Introduction 


The purpose of computing is insight, not numbers. 

— R. W. Hamming 

Synthesis and analysis of control systems for complex physical systems are 
difficult tasks. Complex systems operate in large nonlinear regimes; they often 
comprise many components and contain a large number of state variables. Linear 
systems have been well studied, and there are systematic design procedures for lin¬ 
ear control systems. On the other hand, nonlinear systems are far less understood. 
Nonlinear control synthesis is limited by the lack of available control schemes and 
analysis techniques; traditional design and analysis are seldom applicable to these 
systems. 

Powerful computers can help engineers and scientists unveil the intricacies of 
nonlinear phenomena and can revolutionize the design of complex, high-performance 
control systems. The computational exploration of nonlinear systems has drasti¬ 
cally changed the way we think about the world. This opens the possibility of 
simulating and designing physical systems, without resorting to linear approxima¬ 
tions. 

Yet despite advances in computer technology, nonlinear control systems are 
still difficult to design and analyze. The difficulties arise from the lack of design 
methodologies that actively exploit the special nature of nonlinearities, and from 
the lack of high-level computational tools that can harness the available compu¬ 
tational power. Few nonlinear systems admit closed-form solutions. The ability 
to understand and control these systems requires extensive numerical simulations. 
The burden of the task is on human engineers to carefully design experiments, to 
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distill qualitative information from numerical simulation results, and to use the 
information to design appropriate control actions. This manual, numerical ap¬ 
proach is burdensome and error-prone, and severely limits the design space one 
can explore. 

The modern geometric theory of dynamical systems, pioneered by Poincare 
at the beginning of the century, provides a qualitative way to describe the rich 
dynamical behaviors of nonlinear systems. The geometric representation forms 
a viable substrate for computationally analyzing qualitative aspects of dynamics 
and exploring novel control-synthesis strategies. Programs that employ powerful 
symbolic and numerical techniques have already shown promise in the geometric 
analysis of nonlinear dynamical systems [2]. 

However, the phase-space geometry and topology of a dynamical system are 
hard to describe in a way that allows for efficient computational manipulations for 
the purpose of control design. Nonlinear systems can have extremely convoluted 
phase-space geometries; the complexity becomes much worse as the dimensional¬ 
ity increases. Humans can comfortably picture and manipulate two and three- 
dimensional objects with the aid of graphic, geometric modeling techniques. For 
higher-dimensional systems, however, few visualization and manipulation tools ex¬ 
ist. Automatic modeling, analysis, and design tools are necessary to identify, 
extract, and reason about the spatial properties of phase space. 

1.1 Contributions 

This thesis demonstrates that difficult control synthesis tasks can be automated, 
using a computational workbench consisting of programs that actively exploit 
knowledge of nonlinear dynamics and phase space. These programs combine pow¬ 
erful numerical and symbolic computations with spatial reasoning techniques. The 
thesis contributes to the state of the art in artificial intelligence and control systems 
engineering in several ways: 

• The thesis has developed a qualitative representation for complex behav¬ 
iors of dynamical systems in phase space and a design language for com¬ 
putationally expressing and manipulating these behaviors. The qualitative 
representation captures the gross aspects of dynamics in a relational graph 
of phase-space structure and a set of discrete objects called flow pipes —the 
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equivalence classes of behaviors. The design language describes a control de¬ 
sign task in terms of well-defined geometric, combinatorial operations on the 
flow pipes. This language helps formalize aspects of implicit expert reasoning 
of control engineers in solving control design problems. The representation 
and the language are developed independently of the orders of systems, i.e n 
the dimensionality of phase spaces. 

• The thesis has constructed a computational environment, the Control En¬ 
gineer’s Workbench, integrating a suite of programs that automatically 
analyze and design high-performance, global controllers for a large class of 
nonlinear systems. These programs combine powerful techniques from nu¬ 
merical and symbolic computations with novel representation and reason¬ 
ing mechanisms of artificial intelligence. The two major components in the 
Workbench—the Phase Space Navigator and MAPS 1 —work together to vi¬ 
sualize and model the phase-space geometry and topology of a given system. 
They reason about and manipulate the phase-space geometry and topology 
and search for optimal control paths connecting initial state and the desired 
state for the system. The Workbench represents the result of design and 
analysis in a symbolic form manipulable by other programs, and produces 
a high-level summary meaningful to professional engineers. It also presents 
the result in a graphical form. 

• The thesis presents a novel phase-space design methodology for the syn¬ 
thesis of automatic controllers. The design methodology computationally 
exploits dynamical systems’ nonlinearities in terms of phase-space geome¬ 
tries and topology. It designs a prespecified control law—control reference 
trajectories—for a system by synthesizing the desired shape for phase-space 
geometry dictating trajectory flows. It uses the novel technique of flow pipes 
to group infinite numbers of distinct behaviors into a manageable discrete set 
that becomes the basis for establishing reference trajectories, and navigates 
the system along the planned reference trajectories. The phase-space design 
approach requires powerful computational tools that are able to identify, ex¬ 
tract, represent, and manipulate qualitative features of phase space, and is 
embodied in programs comprised in the Workbench. 

1 MAPS stands for Modeler and Analyzer for Phase Spaces. 
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• The thesis has demonstrated the Workbench and the phase-space design 
methodology in an application of great practical interest: the Workbench 
helped design a high-quality controller for a magnetic levitation system—the 
German Transrapid system. The control system synthesized by the Work¬ 
bench stabilizes a maglev vehicle with much larger initial displacements than 
those allowed in a previous linear design for the same system using classical 
linear feedback technique. The simulation shows that our design outperforms 
the linear design by a factor of 20. Professional control engineers consider 
this result important enough to present at the 31st IEEE Conference on 
Decision and Control [56]. 

1.2 Goal: Mechanizing Control Design 

A real-world control system is a complex closed-loop system with extremely rich 
dynamics. Computation and reasoning are pervasive in the design and operation 
of the controller. Sensors collect a large amount of quantitative information. State 
and parameter estimators infer hidden information about the system from the 
sensed data. The system is modeled with a representation appropriate for further 
analysis and design based on available information. The model is then analyzed 
to extract behaviors that are judged significant for the control objective. To meet 
the control objective, a control law is synthesized to change the natural dynamics 
of the system. 

The domain of automatic control brings together issues of sensing, estimation, 
control synthesis, and control execution. The study of their common themes— 
computation and reasoning —serves as a framework for coherently addressing 
these issues and makes it possible to employ advanced computational techniques 
to drastically improve modern control design. This thesis focuses on the control 
synthesis. 

A control engineer goes through the following design steps to synthesize a 
controller for a given physical system: 

1. Modeling: translate the physics and the constraints of the physical system 
into a quantitative mathematical model. 

2. Analysis: analyze the model of the system. 

3. Design: design a controller for the system. 
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4. Verification: verify the control design. Iterate Steps 2-4 if necessary. 

5. Implementation: implement the control design in software or hardware. 

In the above control design procedure, Step 1 models the system, usually with a 
set of differential equations. Step 2 examines the model and analyzes the behaviors 
of the system. Based on the analysis, Step 3 arrives at a control design according 
to the prespecified control constraints. Step 4 ensures that the design meets the 
control specification. The last step implements the control design to control the real 
physical system. The Steps 2, 3, and 4 are often iterated before a reasonable control 
law is synthesized. Except for very few cases in which analytic-form solutions are 
available, computer simulations are the main tools for analyzing the systems and 
for designing and verifying the controllers. 

The task of control design and analysis requires powerful knowledge rep¬ 
resentation and reasoning mechanisms and computer simulation tech¬ 
niques. The mechanization of this task and the exploration of novel 
control design methods are the subjects of this thesis research. 

1.3 Limitations of Traditional Methods 

We are interested in developing a computational design method suitable for dealing 
with the complexities of real-world nonlinear control systems and in mechanizing 
the design process by computer programs. Traditional methods relying on purely 
analytic or numeric techniques or on linear design theory are inadequate for car¬ 
rying out and automating the design process. 

Systematic techniques for control design of linear systems have been well de¬ 
veloped. Controllers for linear systems can be synthesized via methods such as 
pole-placement or Bode diagrams [35], with satisfactory performance. However, 
complex control systems operate in larger regions where system nonlinearities are 
too significant to ignore. These nonlinear control systems do not lend themselves 
easily to linear (small signal) analysis and design. There have been attempts to 
extend linear control techniques to nonlinear control systems: a piecewise-linear 
control design, also called gain scheduling, approximates a nonlinear function with 
several linear pieces, each of which admits a linear control design. For example, 
in aircraft autopilot design, an airplane system model is linearized around several 



hundred operating points selected within the plane’s operating region. Each oper¬ 
ating point is subject to a linear controller design. However, this type of control 
can be very expensive and complex in both design and implementation. For a 
highly nonlinear system operating in a large region, an enormous number of lin¬ 
earizations and controller designs have to be carried out. The complexity and the 
cost of conventional implementations prohibit many practical applications. 

The computational method developed in this thesis designs a nonlinear control 
system in phase space. Unlike gain scheduling, our method explicitly models and 
actively exploits the nonlinearities in terms of trajectory flows and can carry out 
control design automatically. Although phase space is an extremely useful medium 
for exploring novel nonlinear design methods, the phase-space design approach 
requires sophisticated computational techniques and representational mechanisms 
to make the approach computationally feasible. 

Existing control simulation softwares are inadequate for automatically design¬ 
ing highly complex nonlinear systems. Commercially available programs like MAT- 
LAB and SlMULAB [30] rely on numerical simulations. These programs are, at their 
very best, semi-automatic and serve as interactive design aids to human engineers. 
Although they are equipped with elaborate graphic interfaces, these programs 
provide only fragmented, limited capabilities such as integrations and equilibrium 
determination for performing the simulation task; human users need to prepare the 
simulation and to interpret the result. The specialized control toolboxes embedded 
in these programs are “shallow” expert systems; they lack deep domain knowledge 
and do not have mechanisms for computationally representing and manipulating 
a control design. 

1.4 Characteristics of the Approach of This Thesis 

The computational approach presented in this thesis has the following character¬ 
istics: it can analyze and design nonlinear, global control laws; it is autonomous; 
and its development is independent of the order of a system. 

Our approach exploits the richness of nonlinear dynamics and phase-space rep¬ 
resentations. For example, the programs use the flow pipes to explore a variety 
of control paths otherwise difficult to obtain with just a piecewise linearization. 
Our approach designs control trajectories from a global point of view and thereby 
maximizes the effect of control actions to obtain good performance. 
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The programs embody significant knowledge of dynamical systems theory and 
control theory. Because of their autonomous nature, the programs can perform 
certain tasks that human counterparts have difficulty with, for instance, the task 
of visualizing and manipulating complicated high-dimensional geometric objects. 
The machinery for geometric modeling is developed independently of the spatial 
dimensionality. Hence the programs are not restricted to low-order systems, as 
long as the computational complexity allows. 

These programs are mostly useful in synthesizing high-performance control 
systems that do not lend themselves to traditional design and analysis techniques. 
On the other hand, the qualitative analysis could also assist engineers in exploring 
novel control design methods. Aside from the utility of these programs, this thesis 
formulates the task of control design in a computational way; this formulation 
provides a concrete substrate and a rich domain for exploring and developing new 
reasoning and representation techniques otherwise hard to grasp in traditional 
computer science and artificial intelligence research. 

1.5 What the Control Engineer’s Workbench Does 

Given a model of a physical system and a control objective, the Control Engineer’s 
Workbench analyzes the system and designs a control law achieving the control 
objective. A user typically interacts with the Workbench in the following way. 

The user first tells the Workbench about the system: he inputs a system model 
in terms of a differential equation, parameter values, and bounds on state variables 
for analysis in the form of a phase-space region. The user also tells the Workbench 
about the requirements on the control design: he specifies the desired state for the 
system to settle in, the initial states of the system, the allowable control parameter 
values, and the constraints on the control responses. 

The user then asks the Workbench to analyze the system within the parameter 
ranges of the model. The Workbench visualizes the totality of the behaviors of 
the system over the parameter ranges; it represents the qualitative aspects of the 
system in a data structure and reports to the user a high-level, symbolic summary 
of the system behaviors and, if necessary, a graphic visualization of the phase-space 
qualitative features. 

Next, the user instructs the Workbench to synthesize a control law for the sys¬ 
tem, subject to the specified design requirements. The Workbench searches for 
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the global control paths that connect the initial states of the system and the de¬ 
sired state, using the qualitative description about the system. More specifically, 
the search is conducted in a collection of discrete entities representing trajectory 
flows in phase space. After the global control paths are established, the Work¬ 
bench determines the controllable region of the system and the switching surfaces 
where control parameters should change values. A synthesized control reference 
trajectory consists of a sequence of trajectory segments, each of which is under 
a constant control. For a point-to-point control design, the Workbench further 
constructs a smoothed trajectory connecting an initial state and the desired state. 

The Workbench has designed a high-performance nonlinear control law for an 
electro-mechanical system of practical interest—the stabilization of the German 
Transrapid magnetic levitation transportation system. The details of this design 
will be described in Chapter 5. Here, we briefly illustrate what the Workbench 
does in this application. 

The train and guideway system of the Transrapid is described by a third-order 
nonlinear differential equation: 



Ensuring smooth ride-quality for passengers is one of the top priorities in the 
design of the Transrapid system. Unfortunately, the system is unstable at the 
operating equilibrium for parameter V t at 140 volts. The objective of the control 
design is to synthesize a control law that stabilizes the system at the operating 
equilibrium. 

However, the stabilization control of the train, typically traveling at about 350 
miles per hour along the guideway subject to various disturbances, is a very difficult 
task, for the system operates in a large operating region where nonlinearities come 
from the geometries of the train and the guideway and from the magnetic forces. 
The system requires a high-quality global stabilization design. 

The Workbench is given the system equation in a symbolic form, the ranges 
of parameter values, a bounding box in phase space, and the control objective. 
It analyzes the behaviors of the system for the given parameter values and syn¬ 
thesizes a controllable region about the previously unstable operating equilibrium; 
the region is maximum with respect to the allowable control authority. 
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<equilibrium-points: 

equilibrium 1. (saddle at (140. 0. 200.))> 

<flow-pipes: 

flow-pipe 1. (from *infinity* to *infinity*) 
boundary: (2D-manifold) 
flow-pipe 2. (from ^infinity* to ^infinity*) 
boundary: (2D-manifold)> 

Ccontrollable-region: 
interior: 831 tetrahedra 
boundary: 526 triangles 
volume: 0.848 

range-of-control: (140 300) 
number-of-switching-surfaces: 1 
maximum-displacement-z: 0.00455> 


( c ) 

Figure 1.1: The high-level description of the control design for a maglev system, 
automatically generated by the Workbench: (a) a projection of the controllable 
region showing several synthesized reference trajectories and the operating equilib¬ 
rium; (b) a graphical rendering of the controllable region; (c) a symbolic summary 

characterizing the system and the control design. 
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The Workbench reports the results of the analysis and the control design in a 
graphic and symbolic description in Figure 1.1. The description contains graphic 
displays of the controllable region showing the region boundary, the equilibrium, 
and synthesized control reference trajectories from sampled initial states; it also 
characterizes the behaviors of the system and the controllable region in a symbolic 
form. Although the pictures shown here are in black-and-white, the actual screen- 
display for the region produced by the Workbench is in color. The result of this 
application has been presented to professional control engineers and been favorably 
compared with a previous manual, linear design for the system [56]. 

1.6 Organization of the Thesis 

The rest of the thesis is organized as follows. 

Chapter 2 first gives an overview of the Control Engineer’s Workbench. The 
bulk of the thesis is then divided into three parts: Chapter 3 develops an analysis 
method, MAPS, for constructing qualitative representations of dynamics in terms 
of asymptotic behaviors and equivalence classes, z.c., the flow pipes; Chapter 4 
presents Phase Space Navigator for the synthesis of control laws in phase space. 
The Phase Space Navigator searches for optimal control trajectories in a collec¬ 
tion of flow pipes forming a flow-pipe graph. Chapter 5 applies the machineries 
developed in Chapters 3 and 4 to an engineering problem of practical interest: the 
design of a global, nonlinear controller for the stabilization of a magnetic levita¬ 
tion transportation system. Finally, Chapter 6 describes other related work in the 
computational analysis and design for control systems, and Chapter 7 concludes 
the thesis with a summary and future work in the area. 

Although the domain this thesis studies is control design, readers with knowl¬ 
edge of linear algebra, elementary differential equations, and an introductory to 
linear systems or control will be able to go through the thesis with little difficulty. 
I have included necessary background materials on dynamical systems theory and 
control theory to make the thesis self-contained. 
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Chapter 2 

The Control Engineer’s 
Workbench 


This chapter gives an overview of the Control Engineer’s Workbench—an imple¬ 
mented computational environment for analyzing and designing dynamical control 
systems in phase space. The Workbench, consisting of a suite of programs, allows 
engineers to visualize the phase-space structures of the systems, to manipulate the 
dynamics, and to synthesize controllers. 


2.1 Reasoning About Control Design 

How does a control expert reason about a control design task? A professional 
control engineer uses a body of specialized knowledge to carry out the design. 
The engineer develops insight through an analysis of the physical system, uses the 
insight to explore the design space constrained by control requirements, and makes 
engineering judgment about design choices and trade-offs. A particularly intuitive 
design method is the phase-plane method for analyzing and designing a 2nd-order 
nonlinear control system in a phase plane [35]. 

Let us consider an example of the deformation of an elastic column. Putting 
a weight on top of the column, the column oscillates around its principal axis. 
Putting a heavier weight, the column oscillates and settles to a buckled state 
after a short while. How would an engineer analyze the behavior of the system 
and design a control law to stabilize the buckling motion with the phase-plane 
method? 

The engineer first describes the system with a mathematical model. Compiling 
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the dynamics into a suitable model is a nontrivial task. Discussion of model build¬ 
ing is outside the scope of this thesis. We assume that the engineer is given a simple 
model for the column in the form of an equation of motion, a 2nd-order nonlinear 
ordinary differential equation. Starting with the model and a set of typical pa¬ 
rameter values, the engineer determines all possible responses of the system within 
certain ranges of amplitudes and input signals. The responses with different initial 
conditions from the simulation are then plotted in a two-dimensional plane—the 
phase plane (Figure 2.1). A state, a pair of position and velocity, is plotted as a 
point in the plane. A trajectory consists of a sequence of consecutive states and is 
drawn as a one-dimensional curve. With this graphical plot, the engineer identifies 
the natures and locations of critical points. He then inspects the gross features of 
trajectories: the trajectories approach two critical points in spirals, respectively. 
The engineer groups the trajectories into two classes, each of which corresponds 
to those approaching the same critical point. The phase plane gives a qualitative 
picture of the system responses. The shaded region in the picture is the collec¬ 
tion of trajectories going to the left-hand critical point, one of the buckled states. 
The critical point at the middle of the picture is the unbent state of the column. 
In order for the engineer to interpret the responses, he needs to know things like 
system states, trajectories, critical points, and qualitative regions. 

The buckling motion that leads the column to a buckled state traces out a 
trajectory in the phase plane, for instance, a trajectory that starts at a location in 
the shaded region and evolves to the left-hand critical point in Figure 2.1. Suppose 
the engineer wants the column to return to its original unbent state after the oscil¬ 
lation dies down, instead of letting it get buckled. He introduces a control input, 
a force exerting on the column, and observes how the system responses vary with 
the input. With the pictorial description of the responses, the engineer determines 
a solution curve in phase space as a control law that moves the trajectory to the 
middle critical point. 

In short, the engineer goes through the following steps with the phase-plane 
method: 

• simulate the system extensively with different initial conditions 

• plot the behaviors in trajectories on a piece of paper 

• interpret the result with visual inspection 
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Figure 2.1: The phase plane of the buckling column (from Abraham h Shaw 1988). 

• design a control law to obtain desired behavior. 

The phase-plane method is able to analyze the transient and asymptotic behav¬ 
iors of planar systems, whether linear or nonlinear. Kalman had used this method 
to design switching control laws based on the phase-plane plot [26]. His method 
decomposes a phase plane into discrete linear regions, designs linear control law 
for each region, and matches transients on the boundaries of the regions. 

The phase-plane method is a useful tool for analyzing qualitative responses of 
a control system. It is, however, manual, prohibitively expensive, and confined 
to planes. Although the diagrammatical sketch of a phase plane illustrates just 
the qualitative aspects of the system, to obtain a picture like this requires lots of 
human effort in preparing numerical simulations, collecting the numbers, sketching 
the results, analyzing the trajectory plot, and interpreting it in a qualitative picture 
like Figure 2.1. Worse, this method becomes useless in cases when the order of a 
system is greater than two and the nonlinearity results in convoluted phase-space 
geometry. The mechanization of the) task with autonomous computer programs 
would alleviate many of these restrictions. 
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name : 

- 1 - 

bucliling_column 


equation_of jnotion: 

d : 

h/dt = x 2 



[ dx 2 /dt = —p\Xi — 

P2%\ ~ P3%2 + U 

state_variable: 




state_variable: 

%2 



parameter : 

Pi = 

-2.0 


parameter : 

P2 = 

1.0 


parameter : 

P3 = 

0.2 


parameter: 

U — 

0.0 


boundingibox: 

X\ E 

[-3.0,3.0], 

x 2 e [-4.0,4.0] 


Figure 2.2: Engineer Input: the model for a buckling steel column. 


2.2 A Real Scenario with the Workbench 

The Control Engineer’s Workbench automates a significant portion of the control 
engineer’s design task. In the following session, the Workbench autonomously 
analyzes the buckling motion of a steed column under compression and synthesizes 
a control law to stabilize the motion. The interaction between a human user and 
the Workbench is annotated in explanatory italics for ease of reading. 

• The engineer types in the model of the elastic column buckling 
under axial compressive force and asks the Workbench to analyze 
the system for the given parameter values (see Figure 2.2). 

• The Workbench summarizes the qualitative behaviors of the sys¬ 
tem, in terms of equilibria, stability regions, and flow pipes (see 
Figure 2.3). 

The Workbench also displays the phase portrait of the system 
showing the first flow pipe (see Figure 2.f). 

• The engineer instructs the Workbench to design a control law that 
stabilizes the buckling moti\on, in particular, a law that brings the 
column back to the unbuckled state corresponding to the unstable 
equilibrium at the center of Figure 2.4. He specifies the control 
design requirements (Figure 2.5). 
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<equilibrium-points: 
equilibrium 1. (attractor at 
equilibrium 2. (saddle at (0 
equilibrium 3. (attractor at 


(1.41 (f.)) 
0 .)) 

(-1.41j 0 . ) )> 


<traj ectories: 

<boundary-traj ectories: 
trajectory 1. (from ♦infinity* to (0. 0.)) 
trajectory 2. (from *infinity* to (0. 0.))> 
<connecting-trajectories: 
trajectory 3. (from (0. 0.) to (-1*41 0.)) 
trajectory 4. (from (0. 0.) to (1.4l 0.))>> 


<stability-regions: 
stability-region 1. 

attractor at *infinity* 
stability-boundary: () 
connecting-trajectories: () 
stability-region 2. 

attractor at (1.41 0.) 

stability-boundary: (trajectory 2. trajectory 1.) 
connecting-trajectories: (trajectory 4) 
stability-region 3. 

attractor at (-1.41 0.) 

stability-boundary: (trajectory 2. trajectory 1.) 
connecting-trajectories: (trajectory 3)> 


<flow-pipes: 

flow-pipe 1. (from *infinity* to (-1.41 0.)) 

boundary: (trajectory 2. trajectory 1. trajectory 3.) 
flow-pipe 2. (from *infinity* to (1.41 0.)) 

boundary: (trajectory 2. trajectory 1. trajectory 4.)> 


Figure 2.3: Workbench Output; a symbolic summary of the analysis. 
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Figure 2.4: Workbench Output: the phase portrait of the buckling column showing 
equilibrium points, boundaries of stability regions, and a flow pipe leading to the 
left-hand attractor. 


control.type: 

point_to_point 

goal_state: 

(0.0, 0.0) 

initial _state: 

(-1.0, -3.0) 

range_of .control: [-0.2,0.2] 

-i- 









• The Workbench reports a synthesized control reference trajectory 
(see Figure 2.6 for a phase-space illustration and a script of the 
control law). The global pcrtion of the reference trajectory shown 
in the figure consists of four segments, each of which starts at a 
switching state marked as a small circle. The control parameter 
value is held constant for each segment and is denoted by Ul, U2, 
U 3, or UA. The control IfLW for the reference trajectory is rep¬ 
resented as a sequence of tuples: (time, switching-state, control). 
The control law invokes a local controller when the trajectory en¬ 
ters the vicinity of the goaf The Workbench determines that the 
control response time for the global portion of the reference trajec¬ 
tory is less than 7 seconds. 


2.3 Overview of the Workbench 

In the previous session with the Workbench, the Workbench has demonstrated the 
following capabilities: 

• Deciding what behaviors are significant. The Workbench looks for qualitative 
features like equilibrium points, stability regions, and trajectory flows. 


• Describing the behaviors quality 
bench models the stability 


ively in computational terms. The Work- 
regipns and trajectory flows geometrically. 


• Reasoning about a phase space with flow pipes. 


• Performing the control design autonomously. 

These capabilities are supported in the Workbench by (1) qualitative represen¬ 
tation consisting of dimension-independent constructs, (2) hierarchical extraction 
of the behaviors, (3) modeling and manipulation mechanisms for flow pipes, and 
(4) algorithms implementing the geonjietric., combinatorial, and numerical compu¬ 
tations. 


2.3.1 Design requirements 


We motivate the design criteria for the Workbench with the intended task domain 
of control synthesis and with computational considerations. 
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;; The Synthesized Control Law: specifying the time instance, switching 
;; state, and corresponding control Value for each switching: 

((time 0.) (switching-state #(-l -3)) (control .2)) 

((time .284) (switching-state #(-1.82 -2.71)) (control 0.)) 

((time 1.06) (switching-state #(-1.80 2.49)) (control -.2)) 

((time 2.71) (switching-state #(1.35 1.82)) (control 0.)) 

((time 6.76) (switching-state #(-.00^3 -.0692)) (control *local-control*)) 

Figure 2.6: Workbench Output: the control reference trajectory and the control 
law for stabilizing the column. In the blot, the reference trajectory is drawn in solid 
lines, and the flow-pipe boundaries ol| the uncontrolled column in dashed lines. 
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needs to present the result of 
and to encourage the designers 


• The design of the Workbench Ce lls for a concise representation that captures 
essential features of a control system. The representation should be sensible 
to human engineers and manipulable by other programs. The Workbench 

its analysis and design to human designers 
to interact with the design in a direct way. 
This communication requires a high-level, intuitive presentation of the result. 
Other programs in the Workbench need to efficiently access and manipulate 
the representation. Instead of encapsulating everything about the system, 
the representation should contain only information that is useful for the 
control synthesis task. We have thosen a qualitative representation describing 
a system in terms of its phase-space geometric features for this purpose. 


• The Workbench needs modeling algorithms to efficiently construct the rep¬ 
resentation. The algorithms should identify and extract implicit dynami¬ 
cal properties from numerical explorations and summarize the result in the 
qualitative form. In the implementation, the Workbench internally uses a hi¬ 
erarchy of intermediate representations. The qualitative information about 
the system is extracted in a step-by-step fashion, from local descriptions to 
global ones. 


• The Workbench also needs a Reasoning mechanism to efficiently manipu¬ 
late the representation for synthesizing a control law. The program searches 
through the representation to find feasible control trajectories. The Work- 
bench uses a graph mechanism that manipulates a discrete collection of ob¬ 
jects called the flow pipes. 


2.3.2 Anatomy of the Workbench 

The Workbench analyzes and designs control systems in phase space. In the sce¬ 
nario in Section 2.2, the Workbench locates equilibria and models stability regions; 
it groups trajectories, i.e., behaviors, into the flow pipes; it constructs a dimension- 
independent representation for the phase space; it designs control reference tra¬ 
jectories with this representation. In addition, the Workbench communicates with 
the user in high-level terms through an interface for inputting the model and for 
describing the design and analysis and a graphic module for displaying the result. 

Like the phase-plane analysis, the Workbench is able to reason about and ma¬ 
nipulate dynamics in terms of phase-sjpace geometries. It partitions a phase space 
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into discrete regions. For the previous elastic column example, the program au¬ 
tonomously explores the behaviors of the system and generates the high-level pic¬ 
ture of the phase portrait in Figure 2.4. The picture contains essentially the same 
kind of information as that of Figure 2.1. 

But unlike the phase-plane method, the Workbench provides computational 
means for manipulating the dynamics. It decomposes the phase space into sub- 
regions that can even be globally nonlinear. The shaded region in Figure 2.4 is 
approximated as a geometric polygon The program internally represents the criti¬ 
cal points and geometries of regions iri a data structure that allows other programs 
to manipulate, visualize, and communicate with human users. Because the geome¬ 
try and topology of a system’s phase ^pace are modeled with a simplic.ial structure, 
the representation and reasoning mechanisms for this structure are independent of 
the dimensionality of phase space. 

The Workbench serves as an intelligent assistant to control engineers. The 
components of the Workbench are shown in Figure 2.7: 

• the MAPS program for simulation and interpretation 

• the Phase Space Navigator for control synthesis 

• a graphic program for visualizing the design 


• a user interface for communication with the system. 


The component for model building is inot in the Workbench yet. The programs in 
the Workbench implement various algorithms: symbolic differentiation, numerical 
algorithms on differential equations, humerical integrations, modeling of geomet¬ 
ric structures, clustering of equivalence classes, graph algorithms, etc. The infer¬ 
ence mechanism of the Workbench uses these programs to construct a qualitative 
phase-space structure for representing a system, to check for the consistency of the 
structure, and to reason about and manipulate the representation through a graph 
of flow pipes. 


The flow of computation within the Workbench is illustrated in Figure 2.8. 
Given a model of a system, a bounded phase-space region of interest, allowable 
parameter values, and control objecti\es and constraints, the Workbench performs 
stability and trajectory flow analysis for the system in phase space and interprets 


the result in a 
to synthesize 


phase-space graph. Thq Workbench then explores the control space 
a desired control law subject to the design constraints. It reports 

i 
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Figure 2.7: The structure of the Workbench. 

























the control law specifying reference trajectories and performance properties. The 
control design for steering towards an equilibrium is performed in this way: for 
a point-to-point control, the output is a reference trajectory, whose control law 
consists of a sequence of tuples of time, switching state and the corresponding 
control value; if an initial operating region is given, the output is a controllable 
region geometrically represented as a polyhedral structure. The steering towards 
a limit cycle has not been implemented yet, although its implementation would 
be very much the same in Poincare Sections as the case for steering towards an 
equilibrium. 

The programs in the Workbench 4re written in Scheme, a dialect of LISP. The 
Workbench is implemented on HP 70Q series workstations. For moderately difficult 
problems, for example, the buckling steel column example, it takes the Workbench 
on the order of minutes to arrive at a control law. 
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Figure 2.8: The flow of Computation in the Workbench. 











Chapter 3 

Automatic Ph^se-Space 
Modeling and Analysis — MAPS 


1 his chapter develops and demonstrates a computational method for automati¬ 
cally analyzing qualitative behaviors of dynamical systems in phase space. We 
will discuss the extraction and representation of the qualitative behaviors of a dy¬ 
namical system in terms of a qualitative phase-space structure describing spatial 
arrangement of geometric entities in phase space, a geometric construct called flow 
pipe for modeling trajectory flows, ap analysis algorithm and its implementation 
in the MAPS program, and examples illustrating the mechanism of the analysis. 


3.1 Introduction 

Analysis of dynamical systems via phase-space geometric structures plays an in¬ 
creasingly important role in experimenting, interpreting, and controlling complex 
systems [1, 2]. Nonlinear systems usually fall outside the domain of traditional 
analysis methods, such as Fourier analysis for linear systems. However, most of 
the important qualitative behaviors of a nonlinear system can be made explicit in 
phase space with a geometric analysis [39]. 

MAPS is a program for understanding and representing qualitative structures 
of phase spaces. MAPS combines Numerical, symbolic, and geometric compu¬ 
tations with techniques of spatial reckoning. It embodies theoretical knowledge 
about nonlinear dynamical systems and formulates in computational terms a sig¬ 
nificant amount of informal working knowledge of professional control engineers in 
analyzing complex control systems, in particular, in the form of the phase-plane 
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method [35]. We will illustrate our techniques of extracting and representing the 
qualitative phase-space features with two- and three-dimensional systems. The 
techniques presented in this chapter also apply to higher-dimensional dynamical 
systems. 

Complex systems are often nonlinear and high dimensional. Our theoretical 
knowledge about nonlinear dynamical systems is far from complete. Therefore, 
many engineering applications rely on extensive numerical experiments. A numer¬ 
ical simulation typically generates an immense amount of quantitative information 
about a complex system. To interpret the numerical result and to use the in¬ 
formation for engineering designs, it is essential to develop qualitative methods 
that automatically analyze the system, extract the qualitative features, and rep¬ 
resent them in a high-level description sensible to human beings and manipulable 
by other programs. The representatipn for the qualitative behaviors needs to be 
parsimonious and yet capture the essential properties of the dynamical system un¬ 
der study. The task of control design discussed in the next chapter requires that 
the representation facilitate efficient manipulation in synthesizing new dynamical 
behaviors. We have chosen an equivalence-class based qualitative representation 
motivated by the above requirements on the form and use of the representation. 

Our goal is to develop a class of intelligent and autonomous programs that 
understand the behaviors of complex [systems through their phase-space represen¬ 
tations, synthesize control commands, and affect the physical processes. For ex¬ 
ample, a controller that controls the locomotion of an autonomous walking robot 
would monitor the state of the system, analyze the motion, and command the 
motor to achieve a particular walking behavior. When these systems operate in 
nonlinear regimes and are of high order, their complexities often defy human anal¬ 
ysis. In order for the robot to autonomously execute control tasks with superior 
performance, it is important that th<i robot understands the consequences of its 
control action and synthesizes appropriate control commands to affect the state 
of the system. We are particularly interested in automating the control analysis 

systems that do not lend themselves easily 


and synthesis for a class of nonlinear 


to traditional analysis and design techniques. 
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3.2 Qualitative Phase-Space Structures 

We are interested in representing the qualitative behaviors of dynamical systems 
for control analysis and design. One useful qualitative representation of the phase 
space of a dynamical system is in terms of equilibrium points and limit cycles, 
stability regions, trajectory flows exh biting the same qualitative features, and the 
spatial arrangement of these geometric objects. We call this representation the 
qualitative phase-space structure for a dynamical system and give the following 
definition. 

Definition 3.1 Qualitative phase-space structure; The qualitative phase- 
space structure within a phase-space region of interest, for a structurally stable 
dynamical system addressed in this thesis 1 , consists of 

(a) the number, positions, and stability types of equilibrium points and limit cy¬ 
cles 

(b) the geometric structures of stability regions (basins of attraction) associated 
with the attractors 

(c) the equivalence classes of trajectory flows 

(d) the spatial arrangement of the equilibrium points, limit cycles, stability re¬ 
gions, and trajectory flows. 


3.2.1 Equilibria, limit cycles, and stability regions 

We review some of the basic concepts in dynamical systems theory in order to 
describe the qualitative phase-space structures. Let us consider a single pendulum 
perturbed from its downward resting position. It swings around its vertical axis 
with a smaller and smaller amplitude, and eventually settles to the resting position 
due to friction. We call such a resting state a stable equilibrium point. We say that 
an initial state of the pendulum is in the stability region of the stable equilibrium 
point if the pendulum starting from that state eventually settles to the stable 
equilibrium point. In this example, all initial states that the pendulum starts from 
except for the vertically upward one cjre in the stability region. 


! The scope of the class of systems is disejissed in Section 3.5.1 
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Figure 3.1: Equilibrium points: (^) attractor, (b) repellor, and (c.) saddle. 


In general, the equilibrium points of a dynamical system x' = /(aqu), where u 
is a parameter, are the zeros of the yector field f(x,u) : R n —► R n . Structurally 
stable systems [18] can have equilibritum points of three types: attractor , repellor , 
and saddle , whose local behaviors in phase spaces are shown in Figure 3.1. An 
attractor is an equilibrium point that all nearby trajectories approach in forward 
time. In the pendulum example, the downward resting state is an attractor. We 
call the attractor an asymptotically stable equilibrium point. A repellor is the one 
that repels all nearby trajectories. One can think of it as an attractor in reverse 
time. Trajectories approach a saddle in some directions and leave it in the other 
directions. The vertically upward st^te of the pendulum is a saddle. Start the 
pendulum resting at this position. A slight perturbation would move it from the 
position. On the other hand, the pendulum with just the right amount of energy 
could swing to the upward position 4,nd rest there, although it will be extremely 
rare and will take an infinite amountl of time to do so. We call the repellors and 
saddles unstable equilibrium points. The study of equilibria plays a central role 
in dynamical systems theory. The equilibrium points, together with other phase- 
space entities, partition a phase spac.^ into qualitatively different subregions. The 
attractors are the only kind of equilibria that can be observed physically; they 
characterize the asymptotic behavior^ of dynamical systems. The other classes of 
steady-state behaviors addressed in tjiis thesis are limit cycles. 

The important concept of stability! is associated with the stability regions. The 
collection of trajectories approaching an equilibrium point is called the stable tra¬ 
jectories (or stable manifold ) of the pd>int; and the collection of trajectories leaving 
an equilibrium point is called the unstable trajectories (or unstable manifold) of the 





Figure 3.2: A stability region (basin (j)f attraction) for an attractor in a 2nd-order 
system. 

point. A saddle has stable trajectories along some directions and unstable trajecto¬ 
ries along the other directions. The upion of the stable trajectories of an attractor 
is its stability region, often called thjs basin of attraction for the attractor. The 
stability region of an attractor is open. Every trajectory starting in the stability 
region will be attracted to the attractor, by definition. Attractors cannot be on 
the boundary of the region. The region may be either bounded or unbounded; it is 
unbounded if the boundary contains ijto repellors. An example of a stability region 
is shown in Figure 3.2. 

3.2.2 Trajectory flows 

Most nonlinear systems do not hav< 
form. Analysis of such systems often 
trajectories over a period of time. £ 
number of trajectories, an exhaustb 
space of all the trajectories is compu 
sensitive to uncertainties. To make it 


$ solutions that can be expressed in closed 
resorts to computational exploration of their 
Since a phase space is a union of an infinite 
^e search for plausible control paths in the 
ationally too expensive and numerically too 
possible to find good control paths, we bundle 
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Figure 3.3: A trajectory flow. 

together trajectories with similar qualitative features into a trajectory flow. 

The trajectory flows describe the collective shape and direction of the trajecto¬ 
ries. Figure 3.3 shows, for example, a trajectory flow in the phase space consisting 
of all trajectories ending at the same destination marked by the symbol +. 

With the trajectory flow bundle rbpresentation of the phase space, the explo¬ 
ration is constrained in a much smaller space of equivalence classes of trajectories. 
The expensive fine-grain search in phase space, possibly of high dimensions, is 
avoided. Furthermore, the flow bundle is a natural representation with respect 
to the robustness to the effect of noise and uncertainties, compared with individ¬ 
ual trajectories. The effect of noise on an individual trajectory is modeled by 
thickening the trajectory. 

3.3 Automated Qualitative Analysis of Phase-Space Struc¬ 
tures 

MAPS understands qualitatively different regions and extracts and represents geo¬ 
metric shape information about these regions. Given a dynamical system specified 
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as a system of governing equations, 
scription of the phase-space structure 
description will be used as input to o 
example, for synthesizing control laws 


MAPS generates a high-level symbolic de- 
as the result of the analysis. The high-level 
jther programs for further computations, for 
in the next chapter. 


3.3.1 Theoretical characterization of stability regions 

I 

An essential part of the stability analysis for a dynamical system is the determina- 

i 

tion of the boundaries of stability regions—the stability boundaries. Our algorithm 
for determining stability boundaries is based on a crucial theoretical result char¬ 
acterizing the stability boundaries of a fairly large class of dynamical systems 2 . 
Under certain weak conditions to be discussed in Section 3.5.1, the result of Chi- 
ang et al. [10] shows that the stability boundary for an attractor consists of the 
stable trajectories (stable manifold) of equilibrium points and limit cycles whose 
unstable trajectories (unstable manifold) approach the attractor. This allows us to 
numerically determine a collection of trajectory points on the stability boundary 
through calculations of the stable and unstable trajectories 3 . In planar systems, 
the boundaries consist of curve segments obtained from trajectories. In higher 
dimensions, the boundary surfaces are swept out by a set of trajectories on the 
boundary. The following theorem establishes the theoretical basis for the stability 
analysis of the algorithm. 

Theorem 3.1 Characterization oft stability boundaries for attractors [Chi¬ 
an g et al. [10]]: For nonlinear autonomous dynamical systems satisfying certain 
conditions, the stability boundary for an attractor is the union of stable manifolds 
of equilibrium points and limit cycles j, whose unstable manifolds approach the at¬ 
tractor. 


The three conditions that the systems have to satisfy will be discussed later 
in the thesis. The equilibria that can be on the stability boundaries are saddles 
and repellors. A boundary that sepaiates stability regions for different attracting 
sets is called a separatrix. A degenerate boundary is one that does not separate 
stability regions. 


2 Our algorithm does not use Liapunov 
Liapunov function based methods will be 


fU 


dis 


*The trajectory points can be obtained, 


nctions to construct stability regions. The class of 
cussed in Section 6.2. 

or example, from numerical integration algorithms. 
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Parker and Chua have generalized 
aries for attracting limit cycles 


Figure 3.4(a) shows the stability boundary for the attractor labeled by the 
rightmost +. One of the unstable trajectories of the saddle, denoted by 0, ap¬ 
proaches the attractor as t — > oo. Therefore, the stable trajectories labeled “b” 
of the saddle form the stability boundary for the attractor. Figure 3.4(b) shows a 
stability region with degenerate boundary labeled “b” 

Theorem 3.1 to characterize stability bound- 
The corresponding result remains much the 
same as Theorem 3.1 except that the occurrences of “attractor” in Theorem 3.1 are 
substituted with the words “attracting limit cycle”. Figure 3.5 schematically illus¬ 
trates the stability boundary for an attracting limit cycle. The boundary is formed 
by two repelling limit cycles. Parker apid Chua’s generalization is known to be true 
for many important examples. However, it is largely an intuitive generalization; 

Although our techniques for extracting ge- 
attractors are readily applicable to those of 


no rigorous proof has been given yet 
ometries of stability regions for point 


attracting limit cycles, we shall not discuss them in detail here. 

3.3.2 Extracting and representing shapes of stability re¬ 
gions 

Since a phase space of a nonlinear dynamical system often consists of qualita¬ 
tively distinct points and regions, the) phase-space “shape” of the system refers to 
the geometric information about the structures and spatial arrangements of these 
points and regions. A qualitative ana ysis of the dynamical system determines the 
“shape” of the phase space. The geometric information about these regions, for 
example, is extremely useful in analysing stability properties of control designs for 
complex dynamical systems, such as electric power systems and mechanical control 
systems. 


Definition 3.2 The “shape” of a dynamical system; The “shape” of a dy¬ 
namical system refers to the geometric characterization of the qualitatively distinct 
points and regions of the phase-space representation of the system. 


We need to extract the geometric 
the numerical results about the stab 
niously so as to facilitate further con 
is to be used to estimate the volumes 


information about the stability regions from 
ility boundaries, and represent it parsimo- 
jiputations. For example, the representation 
of the stability regions, to reason about the 
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| (b) 

Figure 3.4: Stability boundaries for an attractor: (a) non-degenerate boundary — 
separatrix; (b) degenerate boundary. 
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Figure 3.5: The stability boundary for an attracting limit cycle. 



spatial relations with the stability bq 
of the regions, to extract geometric ir 


undaries, to compute topological properties 
□.formation about trajectory flows, etc. 


Generation of discrete trajectory points 

The boundary of a stability region is numerically approximated by a collection 
of trajectory points. A one-dimensional boundary curve is approximated by a 
sequence of points on boundary trajectories, a two-dimensional boundary surface 
is approximated by a collection of points on boundary trajectories that sweep out 
the surface, and so on. 

The discrete trajectory points are generated with numerical integration method. 
The method of integrating trajectories from initial points is well-know. However, 
it is difficult to generate a set of trajectory points that evenly populate a surface. 
The difficulty arises from deciding w lere the initial points come from and when 
the integration stops. 

Consider for example a two-dimerjsional stability-boundary surface. Since the 
stability boundary is formed by stablej manifolds of saddles, the trajectories sweep¬ 
ing out the surface can be generated by integrating the trajectories backwards from 
initial points near the saddles in the directions of stable eigenvectors. This inte¬ 
gration method also generalizes to higher-dimensional boundary surfaces. 

In order to maintain relatively even spacings between surface trajectories, the 
program must select a subset from all the possible initial points and terminate an 
integration of a trajectory when the trajectory gets too close to other trajectories. 
We use the following criteria for this purpose: (1) the distribution of trajectory 
points, (2) inter-trajectory distances, and (3) local curvature of the surface. 

Let’s take the surface of a sphere as a stability boundary, and assume a saddle 
is at the north pole and a source at thq south pole. As we integrate backwards from 
the saddle along a set of surface trajectories, these backward trajectories diverge as 
they get near the equator. If the distance between two adjacent trajectories is too 
large, an additional trajectory is needqd to cover the gap. The program backtracks 
in this case: it picks an initial point nekr the north pole, be., the saddle, to generate 
this new trajectory. The local curvatjire of the surface also helps determine if an 
additional trajectory is necessary. As the backward trajectories begin to converge 
at the south pole, some of the trajectories can be terminated for the purpose of 
approximation. 
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A simplicial representation 

Given a set of trajectory points on the stability boundary, a structure on the points 
is needed to make explicit the proxinjiity relationships among these points on the 
boundary. We define the structure on the trajectory points as a graph: the vertices 
are the given trajectory points; the edges join those points that are related to each 
other, where the notion of relation is defined by some metric. Among the different 
structures on the points, a minimal representation—one with fewest edges and 
preserving topological structure—is a polyhedral structure having those boundary 
points as vertices. Furthermore, the resulting polyhedron is contained in the convex 
hull of the boundary points. This polyhedron is not unique. We choose the one 
computed from the Delaunay trianguiation to be discussed shortly. 

A polyhedral structure is a consequence of the so-called simplicial represen¬ 
tation. The elements of a simplicial representation consist of simplices [33]. An 
n-dimensional simplex, or an n-simplex for short, determined by n -f 1 geomet¬ 
rically independent points 4 is the cofivex hull of the points. These points act as 
vertices of the n-simplex. A familiar example of a 0-simplex is a point, a 1-simplex 
is a closed line segment, a 2-simplex is a closed triangle, a 3-simplex is a solid 
tetrahedron, and so on. 

The simplices are basic building blocks from which we construct complicated 
geometric spaces. There are certain ways simplices can be joined to form a new 
object; for instance, two triangles can be glued together at a vertex, or along 
an edge to form a quadrilateral. Given a set of points in a space, a simplicial 
tessellation, or a trianguiation, is a collection of simplices that cover the convex hull 
of the given points and have these poifits as vertices. For example, a trianguiation 
on a set of points on a two-dimensional surface divides up the surface into triangular 
regions. 

Among possible triangulations ovey these points, the Delaunay trianguiation [40] 
offers certain advantages over the others with respect to the quality of triangles 
measured by the regularity of the triangular mesh. The Delaunay trianguiation has 
the property that the circumsphere of any simplex does not contain any vertices in 
its interior. This property makes the Delaunay trianguiation very attractive: the 

A set of points |a?o, ®i, x n } of R^ (N > n) is said to be geometrically independent, if 
and only if the vectors x\ - x 0 , ..., x n - x 0 are linearly independent. Thus two distinct points 
in R form a geometrically independent set; so do three non-collinear points in R? and four 
non-coplanar points in R 3 . 
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Figure 3.6: A Delaunay triangijlation over a set of points in a plane. 

simplices in the triangulation are modt equiangular over all triangulations; in two 
dimensions, the minimum interior angles of its triangles is maximized. Figure 3.6 
shows the Delaunay triangulation on |a set of points. 

We are concerned with deriving a, discrete representation for the geometric 
structure of a phase space that computer programs can efficiently reason about. 
The representation needs to preserve the gross features of shape of the struc¬ 
ture with minimum complexity. Yip has used a point-set representation for the 
phase-space geometries and topologies in analyzing two-dimensional Hamiltonian 
maps [52]. While the point-set representation has been successfully used to recog¬ 
nize patterns of the phase portraits, it is not suitable for the purpose of cutting and 
pasting trajectory flows and synthesiz ng desired phase-space geometry and topol¬ 
ogy* The simplicial structure, on the other hand, is a much coarser representation 
for describing qualitative features of geometric pieces. The simplicial represen¬ 
tation is independent of the dimensions, permits recursive decomposition, and is 
characterized by certain algebraic signatures such as fundamental groups [33]. A 
point-set representation or trajectory representation would be too fine-grain for 
manipulating the shape of phase spaces and would incur unnecessary complexity 
in searching for control paths in phase space. 




Extraction of geometric structures 

The geometric information about a stability region is represented by the polyhe¬ 
dron tightly stretched over the trajectory points on the stability boundary. Ex¬ 
traction of this polyhedral approximation proceeds in two steps: computing a 
triangulation of the convex hull containing the polyhedron, and eliminating exte¬ 
rior triangles. The triangulation is computed with the Delaunay method on the 
set of points which results in a tessellation of the convex hull of these points with 
simplices. The polyhedral approximation is then extracted from the triangulation 
by a sculpture method used in visual information representation [7], followed by a 
more expensive centroid method. Simplices exterior to the polyhedron are elimi¬ 
nated first by heuristic rules and then by testing their membership in the stability 
region with trajectories emanating from the centroids of the simplices. 

We have implemented the Delaunay triangulation in n dimensions that supports 
incremental insertion and deletion of points, adapted from an algorithm proposed 
for planar triangulation [37]. The tim<* complexity of the algorithm grows exponen¬ 
tially with the dimensions, so does the number of simplices in the triangulation. 
For N points in n spatial dimension3, the time and space complexities scale as 
N°( n ) [40]. 

Under the condition that the distribution of boundary points is reasonably 
dense and uniform on the boundary, the polyhedron approximating the stability 
region is contained in the triangulation of the convex hull. Furthermore, the union 
of all the circumspheres of interior simplices is a good approximation for the region. 
In order to extract the polyhedron, simplices exterior to the polyhedron have to 
be eliminated. 


For geometric structures with relatively smooth surfaces, z.e., small curvatures, 
the exterior simplices are eliminated with a relatively fast and inexpensive method 
using the circumsphere heuristics [7]. We observe that only certain type of sim¬ 
plices that satisfy the so-called “visibility conditions” can be exterior ones and are 
therefore candidates for elimination. We associate a value V with each candidate. 
The value V is a goodness measure cf the polyhedral approximation to the true 
boundary contributed by this simple^. The minimum of all U’s of the interior 
simplices is denoted by U mm , a measuje for the overall goodness of the approxima¬ 
tion. We focus the search on candidates satisfying the “visibility conditions” and 
eliminate the exterior simplices by deleting the candidate simplex with the least 
V among all the candidates, until the number of points on the boundary of the 



polyhedral approximation is the same as that of the original set of the boundary 
points and the elimination of this car didate could not further increase V„ 

We explore the properties of interior circumspheres versus those of exterior 
ones. When the distribution of boundary points respects the local geometries of 
the boundary they approximate, the envelope of all the circumspheres of interior 
simplices is close to the true boundary. Consider the maximum distance between 
the faces of a simplex that are on the boundary and its circumsphere: the smaller 
the distance is, the better the simplex approximates the boundary. The distance 
for an interior simplex is smaller than that of an exterior one. We use the inverse 
of this distance as the value of V fpr the simplex. Another possible candidate 
for V is the inverse of the circumsphere radius. Figure 3.7(a) illustrates that for 
smooth shapes the exterior triangles tend to be flat and have larger circumcircles, 
compared with the interior ones. 

We reiterate that the condition oh the distribution of the boundary points has 
to be checked with respect to the shjape of a region, to ensure that the circum¬ 
sphere heuristic rule for elimination Ajvorks. For example, more points are needed 
to approximate the boundary of a region with finger-like narrow parts, as shown 
in Figure 3.7(b). 

The sculpture method is stated as follows: 

Algorithm 3.1 Sculpture method [7\: 

1. Construct a set G of all simplices! in the triangulation and a set H of candidate 
simplices satisfying “visibility conditions”; 

2. Order the simplices of H according to their V values; 

3. Consider the simplex s of H havipg the least V. If the number of vertices of the 
polyhedron formed by the set G (s less than that of the original set of points or 
the removal of s increases the V m i n of H , then delete s from G and H; otherwise, 
return G and terminate; 

4. Update neighbors n t of s in G. Add n l to H. Go to Step 2. 

Visibility conditions: 

• In 2D, simplices with one edge and two vertices on the boundaries. 


38 





Figure 3.7: Circumcircle properties of exterior triangles vs. interior ones: (a) a 
region with smooth shape; (b) a region with narrow-parts. 
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• In 3D, simplices with exactly ore face, three edges, and three points on the 
boundary, or those with exactly two faces, five edges and four points on the 
boundary. 


• In n-dimensions, simplices with one, or two, ..., up to n — 1 proper laces on 
the boundary. 

Each simplex .<? in the triangulatipn has a list of neighbors that share one or 
more common faces with s. When s is deleted from the triangulation, its neighbors 
need to be updated. A simplex satisfying the “visibility conditions” is added to 
H once and is deleted from H at nicest once in the elimination process. Let the 
number of simplices in the set G be M. When the set H is implemented as a heap, 
each insertion or deletion of H costs at most O(log M). The worst-case complexity 


of the algorithm therefore is 0(M log 


M). 


For systems whose stability regiops have interleaved geometric structures or 
holes, the geometric extraction algorithm using the circumsphere heuristic de¬ 
scribed earlier often terminates before! all the exterior simplices are eliminated. In 
this case, MAPS resorts to a more expensive method, the centroid method, that 
determines the membership of each siimplex in the stability region: a simplex is 
classified as being in a stability regiofi of an attractor if the trajectory emanating 
from the centroid of the simplex approaches the attractor in the limit or enters 
another simplex already in the stability region. For example, the buckling column 
discussed in Chapter 2 exhibits a bajided phase space (see Figures 2.1 and 2.4). 
The circumsphere heuristic method terminates before all the boundary trajectory 
points are on the boundary of the polyhedral approximation. The centroid method 
is required to eliminate the remaining triangles. 

The centroid method is implemented by an unraveling algorithm, adapted from 
the cell-to-cell algorithm [22], which links together simplices according to the map¬ 
ping of their centroids under trajectory flows. 


Algorithm 3.2 Centroid method: 

1. Let a be the attractor under consid 
the interior or on the boundaries, 
except for those of A, and C be a 


eration, A be the set of simplices containing a in 
h be the set of all simplices in the triangulation 
set initially empty; 
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2. For a simplex s G G , let p be the 
the centroid of s enters. Delete s 
simplices 5 from 5 to A; otherwise, 


next simplex that the trajectory starting at 
[from G. If p G A, add s and all the reachable 
add s to C and make a pointer from p to s; 


3. If G is empty, return A and terminate; otherwise, go to Step 2. 

The time complexity of centroid hiethod scales linearly with the number of 
simplices in the triangulation. The method requires integration over a time interval 
only once for each simplex. 

The sculpture method and the cehtroid method work together to extract the 
geometric shape of the stability regions. The information about the boundary is 
then used in constructing the flow-pipe representation for the phase space. More 
specifically, the boundary informatioijt is used in labeling the simplices with flow 
directions at the boundary. Details of this will be discussed next. 


3.3.3 Modeling trajectory flow pipes 

The flow pipes are further geometric and dynamical characterizations of stability 
regions. They form the foundation fod the search algorithm for synthesizing global 
control paths. 

The polyhedral modeling of stability regions provides a geometric structure 
for reasoning about the stability boundaries. However, it does not provide in¬ 
formation on how the trajectories flotv within the region—the transient behavior 
very important in high-quality control design. The characterization of trajectory 
flows captures the missing information about the underlying vector-field flows; 
such information is essential for efficiently “looking for” good control reference 
trajectories. 

We introduce the technique of flolw pipes for describing the direction and the 
shape of the trajectory flows. The flow pipes form a discrete representation for 
phase space and yet preserve the fine geometries of the trajectories. The flow- 
pipe modeling, together with that of equilibria and stability regions, characterizes 
the transient and asymptotic behaviors of the dynamical systems in terms of the 
qualitative phase-space structures. Figure 3.8(a) shows a portion of the phase 
space for the Lienard equation x" + 0.5;r' -f x 2 = 1, consisting of all trajectories 
ending at the attractor denoted by the symbol -f. A flow pipe groups this collection 


*A simplex x is reachable from s if there 


is a pointer or a chain of pointers from s to x. 
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of trajectories that exhibit the same qualitative feature into an equivalence class, 
as shown in Figure 3.8(b). 


Definition 3.3 flow pipe: A flow pipe models a collection of trajectories exhibit¬ 
ing the same qualitative features in the phase space. It is an equivalence class of 
trajectory paths, each of which can be continuously deformed to another one, de¬ 
limited by the boundaries of stability regions, trajectories connecting equilibria, or 
boundary surfaces dictated by applications 6 . 

Constructing flow pipes 

Computationally, a flow pipe is constructed by aggregating consistent geometric 
pieces with respect to the way that the flow travels from piece to piece. The 
catalogue of flow types on boundaries of the geometric pieces is a set of consistency 
constraints derived from the underlying dynamics of the vector held. 

Given the simplicial tessellation of the phase space constructed by MAPS, each 
proper face / of a simplex s in the tessellation is classified into one of the three 
types with respect to the direction o ' the flow (see Figure 3.9), according to the 
following labeling scheme: 

j 

1. unidirectional flow: 

label in: Trajectories flow into 4 on /; 
label out: Trajectories flow out of s on /. 

2. bidirectional flow: 

label in-and-out: Trajectories flflw into s on some portion of / and flow out 
of ^ on other portion. 

s 

3. zero flow: 

label tangent: The flow does npt intersect /. Hence the proper face / is on 
the boundary of the flow. 

To determine the direction of thi flow on a face of a simplex, the program 
samples N points evenly distributed in the face and computes the components of 

6 Formally, a flow pipe models a homotbpy equivalence class of trajectory paths delimited 
by the boundaries of stability regions, trajectories connecting equilibria, or application-specific, 
boundaries. If the ends of a flow pipe are squeezed to points, the flow pipe is a path homotopy 
equivalence class of trajectories. 
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(a) 



(b) 

Figure 3.8: Grouping trajectories into flow pipes: (a) a collection of trajectories 
with the same qualitative feature; (b) a corresponding flow pipe. 
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Figure 3.9: Classifying faces of simplices with respect to flows. 


flow vectors along the normals of the face at these points. The stability boundaries 
necessarily separate flows and hence a^e on the flow-pipe boundary. Since the true 
stability boundary is only approximated by polyhedral faces, the above sampling 
method could yield incorrect flow labels on the boundary. The boundary infor¬ 


mation obtained in the characterization of the stability regions described earlier 

i 

is used in labeling these zero flow boundary faces. In addition, the trajectories 
that connect saddles on the boundary and the attractor also separate trajectory 
flows. The simplicial tessellation is refined with the introduction of trajectory 
points sampled along these connecting trajectories. 

A proper face of a simplex is called monotonic with respect to both the flow 
and the simplex, if the face is labeled either type 1 or type 3 in the above labeling 
scheme. A simplex is monotonic with respect to the flow if all of its proper faces 
are monotonic. Similarly, a polyhedron is monotonic if all of its proper faces are 
monotonic. Given two simplices agreeing on a non-monotonic, proper face, the 


common non-monotonic face is candeled in the polyhedron formed by the two 
simplices; see Figure 3.10(a). A non-monotonic, polyhedron continues to grow in 
this way until all of its proper faces are monotonic. The monotonic, polyhedra 
are then aggregated to form flow pipes in such a way as to conserve the flow on 
their proper faces; see Figure 3.10(b). The following algorithm groups simplices 
into smallest monotonic polyhedra and clusters the monotonic polyhedra into flow 
pipes. 


Algorithm 3.3 Flow-pipe construction: 
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Figure 3.10: Construction of flow pipes 
polyhedron by canceling common non- 
polyhedra to form a flow pipe accordit 


: (a) grouping simplices to form a monotonic 
monotonic faces; (b) aggregating monotonic. 
:ig to flow directions at the boundaries. 


1. Label proper faces of each simplex with respect to the flow. 


2. Cluster simplices into monotonic bolyhedra 
Cluster simplices into equivalence 
proper faces of the simplices: twc 
face with the same flow label: in- 


3. Cluster monotonic, polyhedra into 
Cluster monotonic polyhedra into 
on the proper faces of the polyhe^r 
agree on a proper face with consis 


classes with the equivalence relation on the 
simplices are equal if they agree on a proper 
dnd-out. 


flow pipes: 

equivalence classes with the equivalence relation 
a: two monotonic polyhedra are equal if they 
tent flow directions: in and out, respectively. 


4. Order the monotonic polyhedra in each class from Step 3 according to the flow 
direction. Return. 


The resulting flow pipes model th 
pipe is formed by the proper faces lab 


je trajectory flows. The boundary of a flow 
eled tangent of simplices in the pipe. 


Classifying shapes of flow pipes 


Flow pipes capture dynamical shapes. A flow pipe ends at an attractor and/or 
starts from a repellor. Saddles can ndt be in a flow pipe except for the boundary 
of the pipe, since flow pipes further decompose the stability regions. In a bounded 
region, each flow pipe either ends at an attractor or leaves the region. 

Flow pipes exhibit topological shapes. Each flow pipe is classified according to 
its relative shape. A flow pipe around a stable limit cycle is a closed pipe. Most 
pipes are open. The open pipes can be relatively straight or highly wound. One, 
two, or more pipes can side by side wind into a spiral, called n-winding pipes. See 
Figure 3.11. We classify the pipes into the following categories: 

1. open pipes: (a) straight pipes; (b) wound pipes. 


2. closed pipes. 

In the earlier Figure 3.8(b), for exa| 
at its end and forms the stability r 
patterns can occur with more pipes oif 


mple, a single pipe winds into a spiral shape 
egion for the attractor. More complicated 
in higher dimensions. 
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Figure 3.11: Topological shapes of flow pipes. 

3.4 The MAPS Analysis Algorithm 
3.4.1 The algorithm 

We present the following algorithm for analyzing, extracting, and representing 
qualitative features of a dynamical system of any order in the phase space. 

Algorithm 3.4 MAPS analysis: 

(1) Identify qualitative behaviors: 

(a) locate equilibrium points/limit cycles and classify their stability types; 

(b) compute stable and unstable trajectories for each saddle/limit cycle; 

(c) identify those saddles/limit cycles whose unstable trajectories approach 
an attractor; 

(d) the stability boundary for the att 'actor is the union of the stable trajectories 
of those saddles/limit cycles identified in Step (c); 

(e) check if consistency rules are violated. If yes, look for missing equilibrium 
points/limit cycles and go to Step (a). Otherwise, go to the next step. 

(2) Extract geometric structures: 
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(a) for each attractor, collect stability 

(b) tessellate the convex hull of the 

(c) extract a polyhedral approximati(J) 


boundary points; 

boundary points with a triangulation; 
n to the stability region. 


(3) Construct flow pipe: 

(a) refine the triangulation of the sta 

(b) label simplices with flow directioij 

(c) aggregate the labeled simplices in 


(4) Summarize qualitative behaviors 

(a) compile the phase-space data stri 

(b) augment the graph with the geor 

(c) report the graph as the output. 


piiity regions; 
s at boundaries; 
to flow pipes. 


and generate a high-level description: 

cture from Step 1 into a relational graph; 
[iketric. structures from Steps 2 and 3; 


ify 


The set of consistency rules spec 
and are used in the algorithm to autcj 


1. The Existence Rule : Every stabilt 
phase space with multiple attracto 


2. The Separation Rule: Separatricqs 
bounded on all ends. 


The first rule states the existenc 
with multiple attractors. The second 
multiple stability regions. The separat 
two stability regions. 

The first step of our algorithm is 
Parker and Chua [38] for numerical 
nar systems. We have augmented thii 
they suggested to automate the loca 
method [41] used in finding equilibriui|n 
Chua method uses a grid to set up in 
and unstable equilibrium points under 
quire that the initial guesses for sadcil 
seek to automate saddle locating by fo 
most likely places using partial bounc 
fining the initial guesses for the Newto 


the conditions for the stability boundaries 
matically locate missing saddles. 


y region of an attractor has a boundary in a 
rs: 


either form a closed surface or become un- 


e of stability boundaries in a phase space 
rule describes the separation property of 
rices are stability boundaries that separate 


based on a numerical method proposed by 
y determining stability boundaries of pla- 
iir method with the set of consistency rules 
ing of saddles. Since the Newton-Raphson 
points requires an initial guess, the Parker- 
i|tial guesses and is able to find all the stable 
normal circumstances. However, they re- 
es be provided manually by the user. We 
cusing the search for missing saddles on the 
ary information already obtained, or by re- 
n-Raphson method. We want to emphasize 
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Figure 3.12: Th^ flow chart of MAPS. 
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that our algorithm is valid for higher dimensional systems as well and generates a 
symbolic description of the phase-space structure. The Parker-Chua method is de¬ 
signed for numerically analyzing planir systems only. We note that our algorithm 
also finds degenerate boundaries if there are any. 

In the second step, the algorithm extracts the geometric structures of stabil¬ 
ity regions. It uses a triangulation method to tessellate the phase space of the 
system and computes polyhedral approximations to the regions. MAPS internally 
represents a trajectory with a sequence of integration points. It uses an adaptive 
4th-order Runge-Kutta method [41] to integrate trajectories and to generate the 
boundary points that approximate a stability boundary. 

The third step refines the triangi lation from the second step and constructs 
flow pipes, using the information abaut the stability boundaries. The unstable 
trajectories of saddles connect the saddles with the attractors. We call these tra¬ 
jectories the connecting trajectories. The triangulation from Step 2 is refined by 
the insertion of trajectory points from the connecting trajectories. The algorithm 
for the flow pipes clusters simplices into equivalence classes as described in Sec¬ 
tion 3.3.3. 

The last step represents the data structure of the phase-space geometries and 
topologies with a relational graph. The details of each step will be illustrated with 
an example in the following section. 

The flow chart of MAPS is show i in Figure 3.12. The input to MAPS is a 
system of governing equations for a dynamical system. We could also start with 
a set of measured states from experiments. These numerical states could then 
be analyzed and clustered to form phase portraits. The phase space could also 
be reconstructed from a measured time series of just one state variable, using a 
technique called delay coordinates introduced in [36, 48]. 


3.4.2 The main illustration 


We illustrate how MAPS computes the: high-level description of a dynamical system 
with an example. Consider a 2nd-order nonlinear system 


x' = —Sx -f 4x 2 — xy /2 — x 3 
y' = -2 Ay -f xy + u 


(3.1) 


where u is a parameter. For the parameter value u = 0.2, the vector field within 
the region -1.0 < x < 4.0 and -1.0 < y < 4.0 is shown in Figure 3.13. MAPS 
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Figure 3.13: The main illustration: thj? vector field of a 2nd-order nonlinear system. 


the 


analyzes the qualitative behavior of 
actually performs the analysis in a ^ 
trajectories that enter the bounding 
implementation magnifies the origin^ 
side.) The equilibrium points of the 
on f(x,u )—the Newton-Raphson met 
within the region and classifies their 
Jacobians at the equilibrium points: 
and two saddles at (1.05,0.19) and (3 
that the attractors are represented w 
sented with the symbol 0 .). The stab 
then computed by integrating the syst 
in the directions of the stable and un 
respectively. 

Since one of the unstable trajector 
(2.0, 2.0), the stability boundary of the 
of both saddles. Similarly, the stabilit 


system within this region. The program 
lightly larger bounding box to account for 
box again after a short exit. (The current 
bounding box by a factor of 1.5 on each 
system are found by a zero-finding method 
hod. MAPS locates four equilibrium points 
stabilities by inspecting the eigenvalues of 
two attractors at (0.0,0.0952) and (2.0,2.0) 
95, —0.21), all shown in Figure 3.14(a) (Note 
ith the symbol + and the saddles are repre- 
e and unstable trajectories of the saddle are 
em from a small neighborhood of the saddle 
stable eigenvectors backwards and forwards, 


ies of each saddle goes to the attractor at 
attractor consists of the stable trajectories 
y boundary of the attractor at (0.0,0.0952) 



consists of the stable trajectories of ths saddle at (1.05, 0.19), one of whose unstable 
trajectories goes to the attractor. However, within the region of interest there are 
trajectories that leave the bounding box. These trajectories can be conveniently 
thought of as the stable trajectories cf an attractor at infinity. Therefore, the sta¬ 
ble trajectories of the saddle at (3.05,-0.21) form the stability boundary for the 
attractor at infinity, for one of the unstable trajectories of the saddle leaves the 
bounding box. Consistency rules are checked and satisfied. At the end of this step, 
MAPS finds three qualitatively different regions associated with the three attrac¬ 
tors and internally represents the phase-space structure in a data structure: the 
attractors are connected with each other via saddles and associated with stability 
boundaries (Figure 3.14(b)). 

The second step extracts a polyhedral approximation to each stability region 
preserving the gross features of the shape of the region. Consider the stability 
region of the attractor at (2.0, 2.0). The stability boundary is numerically approx¬ 
imated by a collection of trajectory paints (70 in total for this example) relatively 
uniform and dense on the boundary; see Figure 3.15(c). 

A Delaunay triangulation is performed on this set of points. As the result, the 
convex hull of the points is tessellated with triangles; see Figure 3.15(d). Since the 
boundary points are distributed reasonably densely and uniformly on the boundary, 
the polyhedral shape of the stability region is contained in the triangulation of the 
convex hull. MAPS successfully eliminates the triangles exterior to the polyhedron 
with the heuristic method, without resorting to the centroid method. The resulting 
region is shown in Figure 3.16(e). 

MAPS then computes the flow-pipe modeling of the phase space. The triangu¬ 
lation of the stability region obtained above is further refined with additional points 
from connecting trajectories. MAPS labels each triangle in the refined triangula¬ 
tion with flow direction on its faces and clusters these triangles into equivalence 
classes. This results in two flow pipes that form the stability region from the 
second step (Figure 3.16(f)). 

MAPS compiles the data structure for equilibria, stability regions, and flow 
pipes from earlier steps into a relational graph; see Figure 3.17. It reports its 
findings to the user in a symbolic summary in Figure 3.18. 
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(d) 


Figure 3.15: (cont’d) The main illustration—MAPS output: (c) points on the 
stability boundary for one of the attractors; (d) triangulation of the convex hull. 
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(f) 


Figure 3.16: (cont’d) The main illustration—MAPS output: (e) polyhedral ap¬ 
proximation to the stability region computed from the triangulation on boundary 
points; (f) two flow pipes computed from the refined triangulation. The flow pipes 
form the stability region. 
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<equilibrium-points: 
equilibrium 1. (saddle at (3.05 -.21)) 
equilibrium 2. (attractor at (2. 2 .)) 
equilibriums, (saddle at (1.05 .19') 
equilibrium 4. (attractor at (0. .0952))> 

<trajectories: 

<boundary-traj ectories: 
trajectory 1. (from ^infinity* to (1.05 .19)) 
trajectory 2. (from *infinity* to (1.05 .19)) 
trajectory 3. (from *infinity* to (3.05 -.21)) 
trajectory 4. (from *infinity* to (3.05 -.21))> 

Cconnecting-trajectories: 
trajectory 5. (from (1.05 .19) to (0. .0952)) 
trajectory 6. (from (1.05 .19) to (2. 2.)) 
trajectory 7. (from (3.05 -.21) to (2. 2.)) 
trajectory 8. (from (3.05 -.21) to *infinity*)» 

<stability-regions: 
stability-region 1. 

attractor at ^infinity* 

stability-boundary: (trajectory 4. trajectory 3.) 
connecting-trajectories: (trajectory 8) 
stability-region 2. 
attractor at (2. 2.) 

stability-boundary: (trajectory 4. trajectory 3. trajectory 2. trajectory 1.) 
connecting-trajectories: (trajectory 7 trajectory 6) 
stability-region 3. 

attractor at (0. .0952) 

stability-boundary: (trajectory 2. trajectory 1.) 
connecting-trajectories: (trajectory 5)> 

<flow-pipes: 

flow-pipe 1. (from *infinity* to (2. 2.)) 

boundary: (trajectory 3. trajectory 1. trajectory 7. trajectory 6.) 
flow-pipe 2. (from *infinity* to (2. 2.)) 

boundary: (trajectory 4. trajectory 2. trajectory 7. trajectory 6.) 
flow-pipe 3. (from *infinity* to *infinity*) 

boundary: (trajectory 4. trajectory 3. trajectory 8.) 
flow-pipe 4. (from *infinity* to (0. .0952)) 

boundary: (trajectory 2. trajectory 1. trajectory 5.)> 
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Figure 3.18: The main illustratior—MAPS output: a symbolic summary. 



3.4.3 Implementation details 

MAPS is implemented in the Scheme Language. Internally, a point in phase space 
is a vector of the coordinates. An equilibrium point is a list of its position, eigenval¬ 
ues and eigenvectors of Jacobian at tie point, and its stability type. A trajectory 
specifies the type of the trajectory, s ;art and end points, and a sequence of inte¬ 
gration points as a Scheme stream object. The type could be either stable-man 
or unstable-man, standing for stable manifold or unstable manifold; a trajectory 
typed stable-man is integrated backwards, and a trajectory typed unstable-man 
is integrated forwards. A stability region is represented as a list of the attractor, 
the boundary and connecting trajectDries, and the polyhedral approximation. A 
flow pipe is a list of its polyhedral c.pproximation, the boundary, and the start 
and end faces. Table 3.1 summarizes the internal representation of these objects. 
The examples for the stability region and the flow pipe are omitted due to space 
limitations. 

3.4.4 More examples 

Our algorithm also applies to higher-dimensional systems. Consider the following 
3rd-order nonlinear system 

(x' = y 

< y' — z — x 2 — y 

y = 1.0 — ?/ 2 — 0 . 8 T 2 

MAPS locates an attractor at (1.06,0.0,1.12) and a saddle at (-1.06,0.0,1.12) 
within the region -5.0 < z < 5.0, -5.0 < y < 5.0, and -5.0 < * < 5.0, and 
determines that the stable trajectories of the saddle form the stability boundary 
for the attractor. The stability boundary is a two-dimensional surface and is ap¬ 
proximated by a set of relatively evenly spaced trajectories. MAPS then tessellates 
the phase space with tetrahedra and extracts a polyhedral approximation to the 
stability region of the attractor (see Figure 3.19). As described earlier, MAPS uses 
general geometric primitives — the n-simplices — to approximate stability regions 
in n dimensions. 

We have also run MAPS on a system that is not generic and falls outside the 
domain of our method (see discussion in Section 3.5.1). The system is the earlier 
2nd-order example (3.1) in Section 3.4.2 with parameter value u = 0. MAPS 
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Object 

Representation 

Example 

a point 

a vector of coordin 

#(2. 2.) 

ates 

Object 

Representation 

Example 

an equilibrium 
a list: (point, eiger 
(#(2. 2.) (-.< 
(#(1. 0) #(.C 

values/eigenvectors, stability type) 

>5+1.41i -.05-1.41i) 

25 1.)) stable) 

Object 

Representation 

Example 

a trajectory 
a list: (type, start, 

(stable-man *ini 
(#(3.05 -.21) 

end, integration points) 
inity* #(3.05 -.21) 

#[promise 26])) 

Object 

Representation 

Example 

a stability region 
a list: (attractor, b 
polyhedi 

oundary trajectories, connecting trajectories, 
on) 

Object 

Representation 

Example 

a flow pipe 
a list: (polyhedron. 

boundary, start faces, end faces) 


Table 3.1: 


The internal representation of phase-space data objects. 
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Figure 3.19: The analysis of a 3rd-order nonlinear system: (a) projection of stabil¬ 
ity boundary and connecting trajectories in x-z plane; (b) projection of polyhedral 
approximation in x-z plane. 
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locates four equilibrium points withi i the region -1.0 < x < 4.0 and -1.0 < y < 

4.0 and computes the stable and unstable trajectories of saddles. It terminates 
with the following partial description: 

<equilibrium-points: 
equilibrium 1. (saddle at (3. 0.)) 
equilibrium 2. (attractor at (2.1 1.98)) 
equilibrium 3. (saddle at (1. 0.)) 
equilibrium 4. (attractor at (0. 0.))> 

<traj ectories: 

trajectory 1. (from *infinity* to [1. 0.)) 
trajectory 2. (from *infinity* to [1. 0.)) 
trajectory 3. (from (1. 0.) to (3. 0.)) 

trajectory 4. (from *infinity* to (3. 0.)) 

trajectory S. (from (1. 0.) to (0. 0.)) 

trajectory 6. (from (3. 0.) to (2.:. 1.98)) 

trajectory 7. (from (3. 0.) to *in::inity*)> 

<saddle-connection: 

trajectory 3. (from (1. 0.) to (3. 0.))> 

MAPS discovers a saddle connection in the course of determining the relational 
graph of the phase-space structure: the trajectory that serves both as an unstable 
trajectory of the saddle at (1.0, 0.0) and as a stable trajectory of the saddle at 
(3.0, 0.0). The saddle connection is labeled as trajectory 3 in Figure 3.20. MAPS 
concludes that the system is not generic—which also implies structural instability 
for the planar system here—and abandons any further efforts to characterize the 
phase-space structure. Since saddle connections are often precursors to chaos or 
structural instability, we believe that they are important in partially characterizing 
the phase-space structures of chaotic or structurally unstable systems. 

3.4.5 Hierarchical extraction and representation of phase- 
space information 

We have described and demonstrated our algorithm for analyzing a dynamical 
system through successive computations on the system, starting from its system 
equation representation. MAPS generates a high-level description of the dynam¬ 
ical system at the end of the analysis To bridge the large semantic gap between 
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Figure 3.20: A 2nd-order system with a saddle connection. 


the deduced symbolic description and the system equation representation of the in¬ 
put, MAPS has employed multiple intermediate representations for the dynamical 
system, shown in Figure 3.21. 

MAPS extracts the information incrementally, applying a set of operations to 
each intermediate representation. At each level of the representation, implicit 
properties, such as spatial relations, of the system at different scales are made 
explicit and thus can be accessed and manipulated by the operators at the next 
level of the representation. In the orcer of analysis, MAPS first generates a local 
description of equilibrium points, limit cycles, and their eigenstructures. It then 
computes the polyhedral approximations to stability regions and trajectory flows. 
Finally, it pieces together the local information about each limit set to form a 
global picture of the phase-space str icture: the relational graph describing the 
interactions of equilibrium points, limit cycles, stability regions, and flow pipes. 

The internal representation of the phase-space description as a relational graph 
captures the qualitative aspects of the phase-space structure. In the relational 
graph, nodes are attractors and arcs denote the relations between their stability 
regions. Each node has information about the attractor it represents, the associ¬ 
ated stability region, trajectory flows, and their polyhedral approximations, and 
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Figure 3.21: A multi-layered representation for a dynamical system used in MAPS. 


the boundary trajectories and boundary equilibrium points and limit cycles. Each 
equilibrium point has information about its position, eigenvalues and eigenvectors, 
and stability type. 


3,5 Discussion 

MAPS analyzes qualitative behaviors of nonlinear dynamical systems using knowl¬ 
edge about stability and trajectory flows from dynamical systems theory. It ex¬ 
tracts the geometric information from the numerical results and represents it with 
a qualitative phase-space structure. In this section, we discuss the class of dynam¬ 
ical systems to which our method applies, further extensions to MAPS, and the 
role of MAPS in assisting control engineers in visualizing behaviors of dynamical 
systems. 


3.5.1 Scope of the analyses 


We have stated that the theoretical basis of our algorithm for locating stability 
boundaries is the result of Chiang et al. [10]. Chiang et al. give a topological and 
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dynamical characterization of the stability boundaries for a class of autonomous 
dynamical systems and show that the stability boundary of an attractor is formed 
by the stable trajectories of equilibrium points and limit cycles on the boundary. 
This class of dynamical systems is defined by three conditions: 

1. Hyperbolicity: All equilibrium points and limit cycles on stability boundaries 
are hyperbolic 7 . 

2. Transversality: The stable and unstable trajectories of equilibrium points 
and limit cycles on stability boundaries are transversal to each other. 

3. Finite limits of boundary trajectories: Every trajectory on the stability 
boundaries approaches an equilibrium point or a limit cycle as t — > oo. 

Conditions (1) and (2) are generic properties satisfied by almost all dynami¬ 
cal systems. The example in Section 3.4.4 that has a saddle connection violates 
condition (2). Structurally stable systems satisfy the generic conditions (1) and 
(2) (see [19, 46] for a definition of generic properties and [18] for discussions on 
structural stability of dynamical systems). Since many systems encountered in en¬ 
gineering applications are structurally stable, we further restrict ourselves to struc¬ 
turally stable systems for the sake of simplicity 8 . Structural stability implies the 
absence of quasi-periodic, orbits. Condition (3) excludes some structurally stable 
systems. Our method applies to a fairly large class of autonomous dynamical sys¬ 
tems as defined by the three conditions. We note that a periodic non-autonomous 
dynamical system can be converted into an autonomous system [18] and can also 
be analyzed by MAPS. 

In its current implementation, MAPS does not handle systems having chaotic 
attractors, in which stability regions and/or stability boundaries exhibit fractal 
structures. It also assumes that the number of equilibria and limit cycles of a 
system is finite. 

7 Consider the Jacobian of the vector field at an equilibrium point. If all the eigenvalues of 
the Jacobian have non-zero real parts, then the equilibrium point is hyperbolic. Similarly, a limit 
cycle is hyperbolic if none of its characteristic multipliers, the eigenvalues of its Poincare map, 
lies on the unit circle in the complex plane. 

8 This is not to say that structurally unstable systems are not interesting. In fact, Hamiltonian 
systems are important in investigating phenomena in which energy dissipation is negligible. 
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3.5,2 Extensions 
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with a bifurcation analysis, similar to| Abelson’s Bifurcation Interpreter [3]. 


3.5.3 The use of MAPS ip visualization 


The qualitative phase-space descriptidn can assist engineers in designing controllers 
for complex systems. 

We envision that a control engineei 
of phase-space structures when varyin; 


r uses MAPS to observe qualitative changes 
£ control parameters. These changes include: 


1. the birth, disappearance, and movement of equilibrium points and limit cy¬ 
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2. the enlargement and shrinking of stability regions 

3. other phase-space qualitative changes. 


A program has been constructed for g 
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In the next chapter, we automate 
Navigator that autonomously synthes 
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tern whose phase space is high dimen 
the desired control properties are imp 
niques. The topological and dynamical 
and trajectory flows forms the basis 


be 


fd) 


raphically displaying three-dimensional geo- 
es of 3D polyhedral structures with different 
^phics rendering tool, an engineer can visu- 
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3.6 Summary of the Chapter 


We have developed a qualitative meth 
structures of nonlinear dynamical syste 
strate the method. MAPS “looks” at 
regions—the stability regions, models 
extracts and represents the qualitative 
edge of dynamical systems theory to 


od for automatically analyzing phase-space 
ms and have constructed MAPS to demon- 
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trajectory flows with equivalence classes, and 
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spaces. It computes a high-level description of a dynamical system through a com¬ 
bination of numerical, combinatorial, and geometric computations and represents 
the phase-space structure with a relational graph. The qualitative representation 
of dynamical systems and the computational formulation of the knowledge of engi¬ 
neering analysis enable us to build computer programs that automatically design 
high-quality controllers for nonlinear systems. In the next chapter we will use 
MAPS to analyze a collection of phase spaces each of which corresponds to fixed 
control parameter values. 
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Chapter 4 

Automatic Phajse-Space Control 
Synthesis — Pljiase Space 
Navigator 


4.1 Introduction 


This chapter develops Phase Space Navigator, an autonomous system for control 


synthesis of nonlinear dynamical systems. The Phase Space Navigator automat¬ 
ically designs a controller for a nonlinear system in phase space. It generates 
control laws by synthesizing the desired phase-space flow “shapes” for the system 
and planning and navigating the system along good control trajectories in phase 
space. 

The Phase Space Navigator relies on the phase-space knowledge of dynamical 
control systems. It employs the MAPS program to extract and represent qualita¬ 
tive phase-space structures characterizing the qualitative aspects of the dynamics. 
The control synthesis utilizes flow pipes to model phase spaces and to search for 
global control paths. We will present the novel idea of phase-space navigation as 
a paradigm for high-performance nonlinear control design. Algorithms for control 
trajectory planning and tracking will be described. The synthesis method will be 
illustrated with an example of synthesizing anti-buckling control laws for a steel 
column. 
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4.2 Automatic Control Synthesis in Phase Space 

The Phase Space Navigator synthesizes a control system from a geometric point 
of view in phase space. The control of a dynamical system is interpreted as the 
“steering” of the system trajectory, emanating from some initial state, to the de¬ 
sired state by a control signal. 

4.2.1 Overview of the Phase Space Navigator 

The Phase Space Navigator consists of a global control path planner, a local trajec¬ 
tory generator, and a reference trajectory follower. The global path planner finds 
optimal paths from an initial state to the goal state in phase space, consisting of a 
sequence of path segments connected at intermediate points where the control pa¬ 
rameter changes. A brute-force, fine-grain search in high-dimensional phase spaces 
would be prohibitively expensive. High-level descriptions of the phase space and 
trajectory flows provide a way to efficiently reason about phase-space structures 
and search for global control paths. The local trajectory generator uses the flow 
information about the phase-space trajectories to produce smoothed trajectories. 
The trajectory follower tracks the planned reference trajectory, reactively corrects 
deviations, and resynthesizes the reference trajectory if the dynamics of the system 
changes significantly. 

4.2.2 Intelligent navigation in phase space 

The control objective for a stabilization problem is to synthesize a control path, 
along which the physical system can be brought to the goal state and made to stay 
there afterwards under control. We tire particularly interested in cases in which 
physical plants to be controlled operate in large regions where nonlinearities of 
the plants cannot be ignored. Whep no global stabilization control laws can be 
found for the desired state with traditional control techniques, composite control 
paths have to be synthesized. Geometrically, this is the case where initial states 
of the systems are far away from the desired state. Phase Space Navigator takes 
advantage of the underlying dynamics of the phase-space flows to plan the tra¬ 
jectory both locally and globally and switches control at carefully planned time 
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Figure 4.2: The Phase Space Navigator 

4.3 The Phase Space Navigator 

We present a general, two-stage architecture for the autonomous synthesis of non¬ 
linear controllers in phase space. The Phase Space Navigator serves as the pro¬ 
totype for the architecture. It has two main modules: the planning module and 
the tracking module. The planning module synthesizes a desired trajectory, the 
reference trajectory, in phase space based on the dynamics of the nominal model. 
The tracking module follows the reference trajectory and reac.tively corrects local 
deviations. The parameters of the nominal model are estimated at runtime and 
the model is updated. When the dynamics changes significantly, the reference tra¬ 
jectory is resynthesized. Figure 4.2 st ows the interplay between the two modules 
in the context of controlling a real physical plant. 

4.3.1 Reference trajectory generation 

The planning module consists of a glob al path planner and a local trajectory gener¬ 
ator. The global path planner finds coarse global paths from the initial state to the 
goal state, utilizing the flow pipes of trajectories. The local trajectory generator fits 
smooth trajectory segments into the g:obal path by slicing out trajectory segments 
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Figure 4.3: Reference trajectory generation 
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Local trajectory generator 

The local trajectory generator produces smoothed trajectories connecting inter¬ 
mediate points. More specifically, trajectory segments are extracted from the flow 
pipes forming the global path. A trajectory segment is extracted from a flow pipe 
and is continuously deformed to connect intermediate points in the path. The 
local trajectory generator also smooths out sharp interconnections of trajectory 
segments to reduce undesired transients, with local linear controllers or local slid¬ 
ing surface insertion [45]. The noise and certain type of parametric uncertainties 
of the synthesized reference trajectory are modeled with a thickened trajectory. To 
guarantee the robustness the thickened trajectory segment is restricted within the 
flow pipe that contains the original reference trajectory segment. However, the 
method discussed here assumes that the model has no parametric uncertainties 
or structural uncertainties that change the topological structures of phase space; 
these kinds of uncertainties include model order uncertainties or unmodeled high- 
frequency dynamics. 

The planning module generates a plan, i.e., the reference trajectory along with 
the control action for each segment of the reference trajectory. 

4.3.2 Reference trajectory tracking 

The tracking module follows the planned trajectory and reac.tively corrects the 
deviation due to uncertainties in the modeling of dynamics, disturbances, etc. The 
tracking is usually local to the region of the phase space containing the trajectory 
segment; see Figure 4.4. The planned reference trajectory provides the global 
feedforward reference term. The senior and estimator measure the actual state, 
whose difference from the reference term is the local feedback correction term. 

Trajectory tracking and reactive correction 

The Phase Space Navigator tracks a synthesized reference trajectory by sensing 
the state of the system, comparing the state with the reference trajectory, and 
computing the current control action. This tracking requires that the full state of 
the system be observable. If the systeitn stays within the tolerance of the reference 
trajectory, the controller simply looks up a control action in the planned trajec¬ 
tory. If the system deviates from the desired trajectory, a local linear controller is 
switched in to force the system back on track. The trajectory correction is local 
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and reactive. It is designed to ensure the robustness of the system in the presence 
of noise in measurements and small changes in system parameters that do not lead 
to bifurcations. 


Model updating and reference trajectory resynthesis 

The tracking module estimates the model parameters through the observables and 
updates the nominal model. When t ie dynamics of the system changes substan¬ 
tially such that the planned global oath is no longer feasible, the Phase Space 


Navigator calls the planning module 


to synthesize a new reference trajectory. 


4.3.3 The autonomous control synthesis algorithms 

We consider a control system x' = f(x,u), where x £ R n is the state variable 
and u £ R m is the control parameter. The initial state of the system is x s and 
the goal state for the system to reach is x g . The set of admissible values for the 

ning algorithm searches for paths from x s to 
nts, in the phase spaces indexed by control 
algorithm follows the synthesized reference 


control parameter u is {ui}. The plan 
x g , subject to certain design constra 
parameter values {uj. The tracking 
path. 


Algorithm 4.1 The trajectory plam\ 
1 


ing algorithm ; 

Generating phase portraits ai|d parsing them into flow pipes: 

Generate a collection of phase portraits indexed by {«,•}. Parse each of the phase 
portraits into flow pipes. 


2. Constructing flow-pipe graph: 

Aggregate flow pipes to form flowi 


pipe graph. Make x s and x g to be nodes in the 
graph. Weigh each edge of the graph with the cost of traversing the edge. 

Testing reachability of goal x g , 

Cluster the flow-pipe graph into flow-pipe path components. If x s and x g are in the 
same flow-pipe path component, go to Step 6. Otherwise, the goal is unreachable 
with the given set of parameter values {n ? }; continue. 

Searching for partial paths frc^m x s to x g : 

Compute the reachable set R Xs from and the reverse reachable set R Xg from x g in 
the flow-pipe graph. Denote the gajp between R Xs and R x to be G =11/?, 


l 'X s 


/?, 
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5. Tuning control parameter: 

Search for a value of the control parameter u outside such that the corre¬ 
sponding phase portrait has flow pipes that, when added to the flow-pipe graph, 
reduce the value of G. If G is zero, i.e., the gap between R Xs and R Xg is bridged, 
collect the flow-pipe paths from x s to x g and go to Step 7. Otherwise, repeat 
Step 5 until there are no more parameter values to search, in which case return 
the collection of partial flow-pipe paths from x s to x q . 




6. Finding flow-pipe paths from 

If looking for the shortest paths 


x , to x 


9 ' 


from x s to x g , run a standard shortest path 
algorithm on the graph. If searching for all paths p from x s to x g subject to the 
constraint C(p) < C, use the following algorithm 2 . 

(a) let P ’global—paths — 0i Ppartial—paths — { (-^s 5 )}• 

(b) for each path p^ £ Ppartial—paths •> ending at X{ : do 

P Vi = {pi * (xi,x)\x one edge reachable from Xi, x ± x g , 

C(pi) + W(xi,x) < C'}; 

Pp t = {Pi * (Xi,x)\x one edge reachable from x t , x = x g , 

C(pi) + W(xi,x)<CY, 

Ppartial—paths — {Ppartial—path,s Pi) U P g • 5 
Pglobal—paths = Pglobal—paths U Pp ... 

(c.) if Ppartial—paths / 0, goto ( 6 b); otherwise, continue. 

(d) if Pglobal—paths = 0, goto Step 4; otherwise, continue. 

(e) order paths in P g iobal-paths according to their costs and return. 


7. Generating trajectory segments: 

Select a representative trajectory 
pipe path from x s to x g and paste 


segment from each flow pipe forming the flow- 
them together at intersections of flow pipes. 


8. Smoothing trajectories: 

Smooth each intersection of trajectory segments through trajectory homotopy de¬ 
formation. Order the smoothed trajectories with costs and return. 

Algorithm 4.2 The trajectory tracking algorithm: 


2 Note that in the discussion C(pi ) is the 
(xi,xj), and is the path composition op 


cost of a path pi, W(xi, x,j) is the weight of the edge 
erator. 
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1. Sensing and estimating the physical system: 

Sense the state x of the physical system, estimate the model parameters, and 
update the nominal model of the system. If x is in the neighborhood of the goal 
x gi g° to Step 5. If the nominal model changes substantially, signal and go to 
Step 4. Otherwise, continue. 

2. Computing control action: 

Based on the synthesized reference trajectory and the observation from Step 1, 
compute control tracking term u\ rack and correction term u correct . If u correct is 
unobtainable with a local linear jcorrection, check alternative paths to the goal 
from the collection of suboptimal paths. If the goal is unreachable, go to Step 4. 
Otherwise, the control action is u\= u track + u correct ; continue. 

3. Generating control action: 

Tune the control parameter to n track and generate u correct with a local linear 
feedback controller. Drive actuators and go to Step 1. 

4. Resynthesizing reference trajectory: 

Call the planning module to resynthesize a reference trajectory and go to Step 2. 

5. Goal stabilizer: 

Stabilize the system at x g with a local linear controller, subject to uncertainties 
and noise. 

Figure 4.5 shows the flow chart of die planning algorithm. The algorithm takes 
as input the system model, allowable control parameter values, design constraints, 
and desired control objectives. It outputs a synthesized reference trajectory in 
the form of a list of tuples: (time, switching point, control parameter value). 
The planning module has been implemented in Scheme. We will describe how 
the planning module synthesizes a control law for stabilizing a buckling column. 
The tracking algorithm described above provides a conceptual framework, whose 
implementation is an immediate goal of future research. 

4.3.4 Discussion of the synthesis algorithms 

On-line and off-line synthesis 

1 he planning module generates smoot led reference trajectories in the phase space. 
The constraints on response time and the availability of computational resources 
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Figure 4.5: The flow chart of the trajectory planning algorithm of the Phase Space 
Navigator. 
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dictate whether the computation of planning is done on-line or off-line. If the 
planning is done off-line, the synthesized plan is compiled into a table. At the run¬ 
time the controller performs a table lookup for a control action. On-line reference 
trajectory synthesis and tracking require substantial computational power. The 
Supercomputer Toolkit [4] provides a,n ideal platform for experimenting with the 
above control synthesis algorithm on real-time control applications. 

In Step 1 of Algorithm 4.1 for planning, the phase-space structure and flow 
descriptions can be either computed once for all the control parameter values or 
generated on demand in the search fof optimal paths. The latter case will be more 
suitable for on-line synthesis. 

Generalizations of the point-to-point planning 

We present the planning algorithm for a control design problem with one initial 
state and one goal state. The algorithm also applies to control systems of “one 
initial state/many goal states”, “many initial states/one goal state”, or “many 
initial states/many goal states”. 

For a “one initial state/one goal state” or “one initial state/many goal states” 
problem, the Dijkstra’s algorithm for the single source shortest paths in a directed 
graph runs in 0(V 2 ) for a graph with V vertices and E edges [nj. 

A “many initial states/one goal state” problem can be converted to the “one 
initial state/many goal states” problem by reversing the directions of all edges; and 
a “many initial states/many goal states” problem is solved by the Floyd-Warshall 
algorithm in 0(V 3 ) [11]. 


Discrete vs. continuous control parameter spaces 

A version of the synthesis algorithm is presented for the control parameter initially 
taking values from a finite discrete set. An example of such finite-valued control 
systems is the satellite attitude controller, which stabilizes the antenna direction 
of the satellite subject to disturbances by turning on or off high thrust jets. Other 
such examples are switching power regulators. 

The method also applies to control systems with continuous, multiple control 
parameter spaces. The program uniformly samples the continuous parameter space 
at discrete points and then applies the algorithm on these control points. A much 
better approach, however, would search for “land-mark points” that delimit dis- 
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tinct behaviors and then partition the parameter space into equivalent subspaces. 
For example, the Bifurcation Interpreter [3] can be employed to search for the 
bifurcation points defining qualitative changes in dynamics and decompose the 
parameter space into topologically equivalent subparts. 

Suboptimal control paths 

The algorithm finds optimal paths in the flow-pipe graph. However, the subop¬ 
timal paths can be useful when the optimal paths are no longer judged feasible. 
They can be stored in the table in addition to the optimal ones. A controller can 
opportunistically choose among available trajectory paths according to the desired 
properties at a control switching point. 

Comparison with dynamic programming 

To synthesize an optimal control path, dynamic programming discretizes a state 
space and conducts a fine-grain search in the discretization. The cost of the ex¬ 
haustive search could be prohibitive for large regions or in higher dimensions. In 
contrast, the Phase Space Navigator searches for control paths in a manageable 
set of flow pipes. It is possible to use dynamic programming within a flow pipe, 
once the global path has been established. 

4.4 An Example: Stabilizing a Buckling Column 

We illustrate the mechanism of the Phase Space Navigator with a control synthesis 
for stabilizing a buckling steel column. The buckling motion of the column has 
been extensively studied in nonlinear dynamics from a theoretical point of view 
and in structural engineering by practicing engineers. The columns are commonly 
used as strengthening elements in structures; for example, flexible space structures 
use columns in large operating regions. Study in nonlinear dynamics shows that a 
slender steel column can exhibit very complicated dynamical patterns under vari¬ 
ous operating conditions. Therefore, it is important to understand the dynamical 
behaviors of the columns and to devise ways of preventing them from failure. We 
will use the Phase Space Navigator to analyze the behaviors of the column and to 
synthesize a nonlinear controller for stabilizing the column. 
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Figure 4.6: Buckling of a thin elastic steel column due to axial end loads. 

4.4.1 The column model 

A thin elastic column, subject to axial end compressive loads, buckles around 
the principal axis, as shown in Figure 4.6. The nonlinearity is introduced by the 
nonlinear geometric stiffness of the column [31]. Stoker [47] gave a simplified model 
for the column subject to axial compressive force and viscous damping 

rnx" -f ex' -f a,\x + a 3 x 3 = 0, (4.1) 

where cq = A + C — 2 P // and a 3 = B + D — P/l 3 . The state x is the characteristic 
measure of the column deflection from the principal axis and x' is the velocity. 
The column has mass m and length 21. The axial load is P, and the coefficient of 
viscous damping is c. The bending stiffness is modeled by a primary hard spring 
with restoring force Ax -f Bx 3 and a secondary hard spring with restoring force 
Cx + Dx 3 . We rewrite equation (4.1) as a system of first-order equations 

{ x\ = x 2 

, w 3 \ (4.2) 

x 2 = ~ « 3X 3 - cx 2 ), 

where x\ represents the deflection and x 2 represents the velocity. 

The system (4.2) describes the buckling motion of the column and represents 
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Figure 4.7: The buckling column: the phase space of a buckling column showing 
the stability boundaries and connecting trajectories. The horizontal axis x is the 
characteristic measure of the column displacement from its principal axis and the 
vertical axis x' is the velocity. 

only a single mode of vibration. For a long and slender column, vibrations are 
observed to occur primarily in the first mode [32]. 

When the axial load P is less than the critical load P cr iticah the column oscillates 
around the principal axis and returns to the vertical state. It has only one stable 
state corresponding to the attractor at the origin of the phase space. Under a 
heavier load, i.e., P > Pcriticah the column buckles to one side or the other. The 
phase space has a saddle at the origin and two attractors symmetrically arranged 
about the saddle; see Figure 4.7. The positions of the two buckled states depend 
on the external load P, the stiffness of the column, and the length /. The larger 
the load is, the farther away the buckled states are from the principal axis. The 
column breaks when the buckling exceeds certain limit. As P surpasses P cr iticai , 
the column undergoes a pitchfork bifurcation of equilibria [18]: the attractor at 
the origin gives birth to a saddle at the origin and two attractors on two sides. 
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4.4.2 Extracting and representing qualitative phase-space 
structure of the buckling column 

MAPS automatically analyzes the column model in the phase space and extracts 
and represents the qualitative phase-space structure [54]. For parameter values 
cii/m — —2.0, a 3 /m = 1.0, and c/m — 0.2 and phase-space region —3.0 < aq < 3.0 
and —4.0 < x 2 < 4.0, the program reports the following findings and represents 
them internally in a relational graph: 

<equilibrium-points: 
equilibrium 1. (attractor at (1.41 0.)) 
equilibrium 2. (saddle at (0. 0.)) 
equilibrium 3. (attractor at (-1.41 0.))> 

<trajectories: 

<boundary-traj ectories: 
trajectory 1. (from ^infinity* to (0. 0.)) 
trajectory 2. (from *infinity* to (0. 0.))> 

<connecting-trajectories: 
trajectory 3. (from (0. 0.) to (-1.41 0.)) 
trajectory 4. (from (0. 0.) to (1.41 0.))>> 

<stability-regions: 
stability-region 1. 

attractor at *infinity* 
stability-boundary: () 
connecting-trajectories: () 
stability-region 2. 

attractor at (1.41 0.) 

stability-boundary: (trajectory 2. trajectory 1.) 
connecting-trajectories: (trajectory 4) 
stability-region 3. 

attractor at (-1.41 0.) 

stability-boundary: (trajectory 2. trajectory 1.) 
connecting-trajectories: (trajectory 3)> 

The program finds two attractors at (1.41,0.0) and (-1.41,0.0) and a saddle at 
the origin. It generates a high-level description of the phase-space geometry: two 
banded stability regions associated with the two attractors, separated by the sta¬ 
ble trajectories of the saddle at the origin. Figure 4.7 shows stability boundaries 
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Figure 4.8: The buckling column: the flow pipe leading to the attractor on the 
left. 

and connecting trajectories of the two stability regions. Based on the geomet¬ 
ric phase-space representation, the phase space is further modeled with two flow 
pipes formed by aggregating geometric pieces, using the flow-pipe construction 
algorithm. The pipe boundaries consist of the separatrices of the two stability 
regions that approach the saddle and of the trajectories that connect equilibria, as 
described by the program in the following: 

<flow-pipes: 

flow-pipe 1. (from *infinity* to (-1.41 0.)) 

boundary: (trajectory 2. trajectory 1. trajectory 3.) 
flow-pipe 2. (from *infinity* to (1.41 0.)) 

boundary: (trajectory 2. trajectory 1. trajectory 4.)> 

Figure 4.8 shows the flow pipe that ends at the left-hand attractor. 

4.4.3 Synthesizing control laws for stabilizing the column 

We want to stabilize the column at its vertical state to support the axial end 
loads and to prevent it from breaking. Under sufficiently heavy load, the buckling 
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motion of the column leads the column to one of the buckled states and, when the 
states are far away from the principal axis, induces the failure of the column. 

We shall focus on global navigation in phase space that synthesizes global 
reference trajectories leading to the desired goal state. Since the design of local 
linear controllers is relatively well understood, we shall not discuss the trade-offs 
among different designs for linear controllers and shall choose, for the purpose of 
demonstration, a simple linear feedback design. Local controllable regions of such 
linear controllers are quantified, given available control strength, and are used in 
constraining the design of global control paths. 

The goal state of the control is the unbuckled state—the saddle at the origin 
of the phase space that does not have a stability region. We want to synthesize 
a non-zero stability region for the goal and maximize the region, 
column is of the form 

x\ — x 2 

X 2 = - a 3X'i ~ CX 2 ) + U, 

where u is the control. In the model, rnu has the same dimension 
The Phase Space Navigator automatically synthesizes a global 
an initial state x s to the goal state x g . We consider two cases for 
* s 3 . 

Control design I: stabilizing the buckling motion 

The column is initially buckling with sufficient velocity. The initial state x s is far 
away from the saddle x g in phase space in this case. The control parameter u takes 
its values from a small range around zero. 

The Phase Space Navigator constructs a flow-pipe graph from phase spaces of 
u = 0 and u =*local-control*, as shown in Figure 4.9. The graph shows the 
case when the initial state x s is in the flow pipe ending at the left attractor as 
shown in Figure 4.8. Denote by U ct i the local controllable region at the goal x g . 
The two flow pipes of u = 0 intersect U ct i at u-i and u 2 respectively. The edges U\x q 

3 We use this example for the purpose of illustrating the Phase Space Navigator. The sys¬ 
tem (4.3) is actually feedback linearizable. A feedback linearization would cancel the nonlinearity 
of the vector field in x 2 direction. In contrast, our design restrains the control to be less than 
10% of the vector-field strength in the first case and less than half of the vector-field strength in 
the second case. 


The controlled 

(4.3) 

as the force, 
trajectory from 
the initial state 
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Figure 4.9: The flow-pipe graph for the buckling column. 

and u^x~ g represent flows within U cti produced by the local controller. The paths 
(infinity* -> x s -> ui -»• left-attractor) and ^infinity* -> u 2 —► right-attractor) 
represent the two flow pipes of the phase space of u = 0. 

The Phase Space Navigator finds a simple path from x s to x g in the graph, 
consisting of edges x s ui and u\x g . Then the program synthesizes an individual 
trajectory connecting ;r s and x g from the flow-pipe path (xJT[, W^)- To construct 
this desired trajectory, the program deforms an uncontrolled trajectory, emanat¬ 
ing from a? 5 , of the flow pipe represented by the edge Tyuf so that the deformed 
trajectory enters U ct i. The region U ctl is first projected backwards through the flow 
pipes to form a goal projection. Figure 4.10 shows the goal projection from U c u —- 
two thin pipes illustrated in thick solid lines. Then the uncontrolled trajectory is 
deformed towards the goal projection. The programmed deformation is designed 
to push the trajectory towards the nearest goal projection in the controllable di¬ 
rection X 2 ‘ Since the control is more effective when the direction of the vector 
field is relatively orthogonal to the controllable direction ;r 2 , switching points are 
inserted to turn off the controller when the angle between the two directions are 
below some threshold. 

Consider for example the case when the column is initially at the state (-1,-3). 
The control design is specified as: 


87 




Figure 4.10: The goal projection and the deformation of the trajectory. 


control_type: 
goal_state: 
initial_state: 
range_of .control: 


point_to_point 
( 0 . 0 , 0 . 0 ) 
(-1.0, -3.0) 
ue [-0.2,0.2] 


The Phase Space Navigator designed a control law that brings the column back 
to the unbent state. The control law is represented as a list of tuples, each of which 
specifies the time, state, and control value for each switching of control: 

((time 0.) (switching-state #(-l -3)) (control .2)) 

((time .284) (switching-state #(-1.82 -2.71)) (control 0.)) 

((time 1.06) (switching-state #(-1.86 2.49)) (control -.2)) 

((time 2.71) (switching-state #(1.35 1.82)) (control 0.)) 

((time 6.76) (switching-state #(-.0023 -.0692)) (control *local-control*)) 

Figure 4.11(a) shows the synthesized reference trajectory originating at ( — 1, —3). 
The circles indicate places where the control of deformation switches. Each seg¬ 
ment of the reference trajectory delimited by the switching points is under interval- 
constant control, as specified by (71, U2 , U 3, U 4, or a local linear control law as 
the trajectory is in the vicinity of the goal. The global portion of the reference 
trajectory is pushed towards the goal projection with small deformation that is less 
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than 10% of the vector field strength at any state or 0.2 in the normalized unit, 
whichever is smaller. The position x and velocity v = x' of the controlled column 
are plotted against the time t in Figure 4.11(b) and Figure 4.12(c), respectively. 
The control u is shown in Figure 4.12(d). 

When the goal is not reachable with the current deformation, the synthesis 
algorithm tries again with increased control strength if possible. The synthesized 
trajectories could be further optimized with variational techniques on the collection 
of trajectories within the flow-pipe segments [6]. 


Control design II: restoring from the buckled state 

The column is initially near the buckled state (-1.41,0.0). The control objective 
here is to pull the column out of the buckled state and to bring it close to the 
unbuckled state. Since the flow pipe containing the initial state does not intersect 
U c ti in the down stream, a different control strategy must be employed. The 
program uses the following strategy to synthesize a control path. 

Assume the initial state of the column is (—1.5, 0.0). The control objective is: 


control_type: 
goal.state: 
initial_state: 
range_of .control: 


point _to_point 
( 0 . 0 , 0 . 0 ) 
(-1.5, 0.0) 
u G [-1.0, 1.0] 


The uncontrolled trajectory from the initial state would spirally approach the 
attractor in the clockwise fashion. The control is exerted in such a way as to swing 
the trajectory away from the buckled state to approach the goal projection of the 
saddle. When the trajectory intersects the goal projection, the control is switched 
off so that the system slides along the uncontrolled trajectory. As soon as the 
system enters U ct i, the linear controller is switched in. The control strength is less 
than half of the vector field strength or 1.0 in the normalized unit, whichever is 
smaller. The synthesized control law is specified as: 

((time 0.) (switching-state #(-1.5 0)) (control -.187)) 

((time .001) (switching-state #(-l.5 .000187)) (control .187)) 

((time 1.65) (switching-state #(-1.24 -.0011)) (control -.289)) 

((time 3.19) (switching-state #(-l.66 .0153)) (control .638)) 

((time 5.92) (switching-state #(-.527 -.00186)) (control -.454)) 

((time 7.74) (switching-state #(-2.02 .0473)) (control 1.)) 
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(a) 



Figure 4.11: The synthesized control law that stabilizes the buckling column: (a) 
the reference trajectory that leads to the unbuckled state corresponding to the 
saddle at the origin. The column is initially buckling with sufficient velocity; (b) 
the position x of the column plotted against time t. (to be continued) 
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((time 8.08) (switching-state #(-1.75 1.47)) (control 0.)) 

((time 11.1) (switching-state #(-.049 .00163)) (control *local-control*)) 

The corresponding reference trajectory is shown in Figure 4.13(a). The control 
parameter changes at places marked by circles. The reference trajectory consists of 
eight trajectory segments labeled by interval-constant control U 1, f/2, ..., (77, and 
a local control law. The first segment ((/1) and the last segment (local control) are 
very short in length and thus are invisible in the figure. The position x and velocity 
v = x' of the controlled column are plotted against the time t in Figure 4.13(b) 
and Figure 4.14(c), respectively. The control u is shown in Figure 4.14(d). 

Control design III: local control near the goal 

A simple local linear controller of the form u local _ linear = -k- l x 1 - k 2 x 2 is used to 
stabilize the system around the goal (see [55] for details). Given the maximum 
control strength 0.2 for the linear controller, we choose an over-conservative region 
for U ct i determined by 

Ucti — {(^ 1 ,^ 2 ) | 10.835a?! + 0.550x21 — 0-04, |0.797xi — 0.605a:2| < 0.04}. (4.4) 
This U ct i is used in both of the above cases. 

4.4.4 The phase-space modeling makes the global naviga¬ 
tion possible 

The qualitative description of the phase-space structure and the geometric mod¬ 
eling of the trajectory flows provide a “map” for navigating system trajectory to 
the goal in phase space. 

The directed graph constructed from the flow pipes is used in searching for 
global paths and in determining whether the goal is reachable. The flow pipes 
also make it possible to characterize the more microscopic spatial and temporal 
relations between the current state and the goal state. For example, the collection 
of simplic.es near the goal forms a neighborhood of the goal that determines whether 
or not the trajectory has missed the goal in the flow pipe. The deformation of 
global path segments is constrained by reasoning about the spatial relation with 
flow pipes. 
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(d) 


Figure 4.14: (cont’d) The synthesized control law for restoring the column from 
the buckled state: the velocity v(= x') of the column and the control signal u for 
controlling the column are plotted against time t in (b) and (c), respectively. 
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4.4.5 Other control problems 

The stabilization for a buckling column is closely related to the pole-balancing 
problem, where the control objective is to balance a pole at the vertically upward 
state. In the phase space of the pole model, the saddle is the upward state of 
the pole and the two attractors on two sides of the saddle correspond to the same 
downward state of the pole. 

The control objective for pole-balancing is to bring the pole to the upward 
state from any other positions, even when the pole is initially hanging downward. 
The problem has been studied as a standard testbed for many control designs, 
most of which focus primarily on linear feedback control in linear regimes. With 
our phase-space navigation strategy, a globally stable reference trajectory can be 
synthesized, similar to the one for the buckling column discussed earlier. 

The optimal design of such trajectories has many practical implications. For 
example, the efficiency of cargo handling work at shipyards depends largely on the 
operation of overhead cranes [44]. Figure 4.15 illustrates such a crane driven by a 
trolley drive motor and a hoist motor. The bottleneck of the ship unloading oper¬ 
ation is the transfer of cargo from a ship to waiting trucks; therefore, minimizing 
this transfer time brings about a large cost saving. The planar motion of the load 
hanging at the moving trolley can be modeled as a swinging pendulum with a mov¬ 
ing point of attachment and a hoisting cable. The trajectories of such overhead 
cranes could be optimally designed with the Phase Space Navigator. This would 
be an important step towards the full automation of cargo handling work without 
a crane operator. 

4.5 Summary of the Chapter 

We have developed an autonomous control synthesis method, the Phase Space 
Navigator, for synthesizing controllers for nonlinear systems in the phase spaces. 
We have discussed primarily finite-valued control systems as examples to illustrate 
the method and have noted that the method also applies to control systems with 
continuous parameter spaces. 

The Phase Space Navigator synthesizes global reference trajectories using knowl¬ 
edge of phase-space structures provided by the modeling and analysis program. A 
phase space is modeled with flow pipes that are collections of trajectories hav¬ 
ing the same qualitative behaviors. The flow pipes provide a way to efficiently 
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Figure 4.1o. The mechanical model for an overhead crane unloading cargo from 
ships, reproduced from [Sakawa&Shindo, 1982]. The crane is equipped with a 
trolley drive motor and a hoist motor. The planar motion of the load is modeled 
as a swinging pendulum hanging at the moving trolley. 

search for and reason about global structures of phase spaces. The global nonlin¬ 
ear control synthesis becomes a graph search problem with this representation of 
the phase space. We have shown how global reference trajectories can be automat¬ 
ically designed with the example of stabilizing a buckling column. With the novel 
idea of grouping an infinite number of trajectories into a manageable collection of 
flow pipes, the difficult task of synthesizing a nonlinear controller is formulated 
as a computational problem that requires a combination of substantial numerical, 
symbolic, and combinatorial computations and spatial reasoning techniques. 

The leverage of the automatic synthesis method comes from applications whose 
operating regions are grossly nonlinear and on which high-performance global con¬ 
trollers are impossible to synthesize with traditional techniques. In the next chap¬ 
ter, we will demonstrate the application of the Phase Space Navigator to the design 
of a maglev controller. Other potential engineering applications include large flex¬ 
ible space structures, robot manipulator planning, satellite attitude control, and 
switching power regulators. 
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Chapter 5 

An Application: Design of a 
Maglev Controller 

This chapter studies a magnetic levitation system and applies the machineries of 
the Control Engineer’s Workbench developed in Chapters 3 and 4 to the design of 
a stabilizing controller for the system. 

The Workbench has helped synthesize a global, nonlinear controller for the 
nominal model of a German maglev system 1 . We describe the systematic state- 
space design method for determining the global switching points of the controller. 
The synthesized control system can stabilize the maglev vehicle with large initial 
displacements from an equilibrium. The simulation shows that our nonlinear con¬ 
troller possesses an operating region more than twenty times larger than that of 
the classical linear feedback design for the same system. 

5.1 Introduction 

Magnetically levitated trains provide a high speed, very low friction alternative to 
conventional trains with steel wheels on steel rails. Several experimental maglev 
systems in Germany and Japan have demonstrated that this mode of transporta¬ 
tion can profitably compete with air travel. More importantly, maglev transporta¬ 
tion can ease traffic congestion and save energy [14, 28, 50]. 

! The maglev stabilization problem was suggested to the author by Professor Richard Thornton 
of Electrical Engineering at MIT, who is leading the MIT Maglev Consortium. The result of the 
work described here appears in a paper [56] to be presented at the 31st IEEE Conference on 
Decision and Control, December 1992. 
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Maglev transportation uses magnetic levitation and electromagnetic propul¬ 
sion to provide contactless vehicle movement. There are two basic types of mag¬ 
netic levitation: electromagnetic suspension (EMS) and electrodynamic suspension 
(EDS). In EMS, the guideway attracts the electromagnets of the vehicle that wraps 
around the guideway. The attracting force suspends the vehicle about one cen¬ 
timeter above the guideway. In contrast, the EDS system uses repulsive force, 
induced by the magnets on the vehicle, to lift the vehicle. 


5.2 The Maglev Model 


An attractive system such as the EMS system is inherently unstable. We consider 
the control design for stabilizing an EMS-mode train traveling on a guideway—a 
simplified model for the German Transrapid experimental system. The Transrapid 
system is schematically shown in Figure 5.1. It uses attractive magnetic forces to 
counterbalance gravitational forces. 

The state equations for the magnetically levitated vehicle and the guideway are 
described by 



dx z( V, — Rx) | xy 


= 9 

= y 


Loz 0 
__ LqZqX 2 
2 mz 2 



where the state variables x , ?/, and z represent coil current in the magnet, vertical 
velocity of the vehicle, and vertical gap between the guideway and the vehicle, 
respectively. The control parameter is the coil input voltage V t . The other param¬ 
eters are the mass of the vehicle m, the coil resistance R , the coil inductance L 0 
and the vertical gap z 0 at the equilibrium, and the gravitational acceleration g. 
Details of the derivation of the model are discussed in [49]. The nonlinearities of 
the system come from the nonlinear inductance due to the geometry of the magnet 
and the inverse square magnetic force law. 

The system has one equilibrium state at which the magnetic force exactly 
counterbalances the force due to gravity and the vehicle has no vertical velocity and 
acceleration. However, the equilibrium is a saddle node which is not stable. The 
control objective, therefore, is to stabilize the vehicle traveling down the guideway 
and maintain a constant distance between the vehicle and the guideway despite 
any roughness in the guideway. The available control input is the coil input voltage 
Vi m the model (5.1). We further assume that V { is produced by a buck converter 
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Electromagnetic Maglev 
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Air Gap 
Approx 3/8 



Force j Force 
Guideway 



LIFT 


( a ) (b) 

Figure 5.1: EMS maglev system for high-speed ground transportation, represent¬ 
ing a simplified drawing of the German Transrapid design, (a) Electromagnetic 
suspension (from [MTAC Report 1989]); (b) Detail of a suspension magnet, super¬ 
imposed on the field distribution (from [Eastham 1989]). 
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capable of delivering any voltage from 0 to 300 volts. 

A linear control design for the maglev system described in [49] uses the pole- 
placement method. The system is first linearized around the equilibrium. The 
linearized system has unstable poles, i.e., the poles in the right-half of .s-plane. 
A linear feedback is introduced to move the poles to the desired locations in the 
left-half of the 6-plane. Such a control design can bring the system back to the 
equilibrium with an initial displacement of up to 0.2mm from the equilibrium. The 
linear controller saturates at the beginning for larger initial displacements. This 
is because the linearized model no longer approximates the original system well 
in regions far away from the equilibrium. A global, nonlinear control law such as 
a bang-bang control that respects the nonlinearity of the system must therefore 
precede the linear feedback control. However, the real challenge for the nonlinear 
design is to determine the global control law specifying, for instance, the switching 
points. 

5.3 State-Space Control Trajectory Design 

We will describe a nonlinear control design—a switching-mode control—in state 
space for the maglev system with large initial displacements from the equilib¬ 
rium. We will show that this controller can be automatically designed with the 
Workbench comprising MAPS and Phase Space Navigator developed earlier. The 
nonlinear controller brings the system to the vicinity of the equilibrium and then 
switches to the linear controller. 

For the purpose of demonstration, we assume that the vehicle is displaced 
from the equilibrium in the direction further away from the guideway. We will 
concentrate on the global design of the control reference trajectories and assume 
that a linear feedback design is available as soon as the system enters the capture 
region of the linear controller. 

5.3.1 Modeling state-space geometry 

The global control law is designed by analyzing and modeling the state-space 
geometry of the system. The Workbench explores the state space of the system 
and characterizes the state space with stability regions and trajectory flow pipes. 
It composes the state spaces for different control parameter values and uses flow 
pipes to synthesize a composite state space. 
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The state variables x , //, and z in the model are scaled by 1, 10 3 , and 2 X 10 4 , 
respectively. The parameters of model are assumed to be: L 0 = O.lh, z 0 = 0.01m, 
R — ID, m = 10000%, and g = 9.8 m/sec 2 , typical of a large vehicle lift magnet. 
Assume the power supply delivers 140 volts, be., V t = 140, at the equilibrium. The 
Workbench explores the state space in a region bounded by the box {(x,y,z)\x 6 
[0,400],?/ 6 [—300,350], z G [0,600]} and finds the following equilibrium point: 

saddle: #(140. 0. 200.) 

eigenvalues: -17.004+22.963i 
-17.004-22.963i 
24.007 

eigenvectors: #(.23604 .97174 0) 

#(.51331 -.55588 .65384) 

#(.30157 .73255 .61027) 

With the information about the stable eigenvectors of the saddle, the Work¬ 
bench computes the stable manifold of the saddle, a two-dimensional surface. The 
Workbench generates a set of trajectories evenly populating the stable manifold to 
approximate the surface. The trajectories are obtained by backward integrations 
from initial points in a small neighborhood of the saddle. This neighborhood lies 
within the plane spanned by the stable eigenvectors of the saddle. 

Figure 5.2 shows the trajectories on the stable and unstable manifolds of 
the saddle in the //^-projection of the state space. The stable manifold is two- 
dimensional and the unstable one is one-dimensional. The stable manifold sep¬ 
al ates the state space into two halves: trajectories in the upper-half approach 
0 —> oo along one of the unstable trajectories, corresponding to the case in which 
the vehicle falls off the rail; and trajectories in the lower-half approach z — 0 plane 
along the other unstable trajectory, corresponding to the case in which the train 
collides with the rail. 

5.3.2 Synthesizing a global stabilization law 

For an initial displacement above or below the equilibrium, the uncontrolled system 
will follow either a trajectory traveling upwards with increasing z and leaving the 
bounding box or a trajectory traveling downwards and hitting the z = 0 plane. To 
stabilize the system at the equilibrium, it is necessary to synthesize a new vector 
field on both sides of the stable manifold so that trajectories travel towards the 
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Figure 5.2: Stable and unstable manifolds of the saddle for V t = 140 (yz- 
projection). 

stable manifold of the saddle in the new vector field. We consider only the top-half 
here. 

The Workbench first explores the state space of the model for different values 
of Vi and concludes that the larger the V t is, the further away the stable manifold 
is from the z — 0 plane. For V{ = v > 140, the region sandwiched by the stable 
manifold of Vi = 140 and that of Vi = v has the desired property—the vector field 
°f Vi = v in this region is pointed towards the stable manifold of V = 140. When 
v = 300, the region is maximized. 

Similarly, the Workbench finds that the model with V t = 300 has a saddle node 
at (300., 0., 428.o7). The stable manifold of the saddle has a similar structure as 
that of Vi = 140 case, but with larger * coordinate. The y ^-projection of the stable 
manifold and unstable trajectories for Vi — 300 is shown in Figure 5.3. Appendix A 
contains the information about the saddle for K = 140 and V { = 300. 

Figure 5.4 shows the ^-projection of the sandwiched region discussed above. 
The region is bounded by three pieces of triangulated surfaces, the top boundary 
shown in Figure 5.5(a), the side one in Figure 5.5(b), and the bottom one in 
Figure 5.5(c), respectively. The top boundary represents the stable manifold of 
the saddle for Vj — 300 and the bottom one for Vj = 140. The trajectories of 
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Figure 5.3: Stable and unstable manifolds of the saddle for V t = 300 (yz- 
projection). 

V; = 300 flow into the region from the side boundary in Figure 5.5(b) and leave 
the region at the bottom boundary in Figure 5.5(c). There is no flow across the 
top boundary in Figure 5.5(a). 

The Workbench determines that the above region consists of two trajectory 
flows, the one for Vi 300 that flows into the region on the side boundary and 
pierces through the bottom boundary, and the one for V t = 140 on the bottom 
boundary that approaches the desired equilibrium in the limit. With the flow-pipe 
characterization of the state-space trajectory flows, the Workbench searches for 
control trajectories in this set of flow pipes and finds a sequence of flow pipes that 
lead to the desired goal: the composite of the trajectory flow for Vi = 300 and the 
flow for Vj 140, glued together at the stable manifold of Vi = 140 represented by 
the bottom boundary of Figure 5.5(c). As a result, all the trajectories of Vi = 300 
within the region can be brought to the equilibrium by switching to Vi = 140 
as soon as the trajectories hit the bottom boundary. We call this region the 
controllable region for the system, as shown in Figure 5.4. 

Since the bottom boundary is an approximation to the true stable manifold, 
the trajectories of Vi = 140 on the boundary can only get close to the desired 
equilibrium. The closeness depends on the quality of the manifold approximation. 
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Figure 5.4: The sandwiched region in ^-projection. 
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Figure 5.5: The boundaries of the sandwiched region in yz-projection: (a) top 
boundary; (b) side boundary; (c.) bottom boundary. 
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As the trajectories enter a small neighborhood of the equilibrium, we use a linear 
feedback controller, such as the one discussed in [49], to stabilize the system at 
the equilibrium. Figure 5.6 shows the synthesized control reference trajectories 
originating at different initial displacements from the equilibrium: 1mm, 4mm, 
4.5mm, and 5mm. The controller is able to recover from the first three initial 
points that are within the region. The last point is outside the region and thus 
uncontrollable; the current in the magnet can not build up fast enough to keep 
up with the ever-increasing airgap. The complete print-out of the synthesized 
reference trajectories showing the switching points is in Appendix B. 

For example, the control law for the initial displacement of 1mm is specified 
as: 


((time 0.) (switching-state #(140. 0. 220.)) (control 300.)) 

((time .0134) (switching-state #(163. 4.44 221.)) (control 140.)) 

((time .127) (switching-state #(139. .0789 202.)) (control *local-control*)) 

This control law specifies the time instance for each switching, the switching 
state, and the corresponding control value during the following time interval. The 
corresponding reference trajectory consists of trajectory segments delimited by 
these switching states. 

5.3.3 Evaluating the control design 

The synthesized global control law is a switching-mode one that changes the control 
parameter at the switching surface—the stable manifold of V- = 140. It is able to 
bring trajectories originating from any states within the controllable region to a 
local neighborhood of the saddle. 

The responses of the controller with respect to the four different initial displace¬ 
ments are shown in Figures 5.7 and 5.8. The vertical axis of each graph represents 
state variables x , ?/, and 2 as in the maglev model (5.1), one for each curve, and the 
horizontal one represents the time. For all the controllable initial displacements, 
the controller is able to bring the system back to the equilibrium with errors less 

than 0.2mm in displacement—a distance within the capture range of the linear 
feedback controller. 

By exploring the state-space geometries of the maglev system, the Workbench 
is able to automatically determine the switching points for the global controller. 
The linear feedback controller can recover from only displacements of less than 
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(b) 


Figure 5.6: The synthesized control reference trajectories originating from four 
different initial states, together with the controllable region: (a) ^-projection; (b) 
zrc-projection. 
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0.2mm. The global controller has significantly enlarged the operating region of the 
linear controller. With the geometric representation of the controllable region in 
state space, the Workbench precisely determined that the maximum recoverable 
displacement is 4.55mm. Our simulation geometrically explained the observation 
in [49] that the vehicle would fall off the rails with a displacement of 5mm or larger. 

Many issues remain to be addressed in order to make the control design prac¬ 
tical. Since our control law is designed with the nominal model of the maglev 
system, the effect of uncertainties in the model and of noise in the environment on 
the design needs to be studied in future research. The design can also be optimized 
with respect to response time. 

5.3.4 Visualizing the design 

The synthesized controllable region is modeled with a polyhedral structure. This 
structure can be presented to engineers in a visual way. J. Choi has implemented 
a graphic rendering program for this purpose. Figure 5.9 shows a picture of the 
graphically rendered controllable region. The graphical presentation facilitates in¬ 
teractive, incremental modifications to the design in terms of the geometry. For 
example, one might want to directly manipulate the controllable region by tuning 
knobs corresponding to parameters. The geometric representation helps the de¬ 
signers to develop intuitions and to visually explore the effects of certain design 
choices. 

5.3.5 Implementation of the controller 

The control design described above has been computationally simulated only. How 
will this design be implemented on the real system and used in real time? 

The control law specifying the switching surfaces in state space can be com¬ 
piled into a table. The control execution will be a table lookup and a geometric 
inequality test. At each step of the execution, the state of the system is sensed 
and checked against the switching surface. If the state is on the switching surface, 
the corresponding control value for the next time interval is read from the table 
and applied to the physical system. The implementation does not have to be very 
different from that for a dynamic programming one. 
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Figure 5.9: The controllable region for the maglev control design, 
the light source straight on. Hidden faces and lines are removed. 


rendered with 
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5.4 Summary of the Chapter 


We have described the automatic design of a high-quality controller for stabilizing 
an EMS-mode maglev vehicle traveling above a guideway with the Workbench. 
The synthesized control reference trajectory consists of a sequence of trajectory 
segments, connected at intermediate points where the control voltage switches. At 
run-time, the controller will cause the system to track the reference trajectory and 
reactively correct local deviation from the desired trajectory. We have illustrated 
the state-space method for designing the global switching points of the nonlinear 
controller. The simulation showed that our design allows the maglev train to oper¬ 
ate with much larger disturbances in the airgap than the classical linear feedback 
design does. 
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Chapter 6 


Related Work 


This chapter surveys relevant work in areas of simulation and analysis for complex 
dynamical systems and nonlinear control synthesis methods. Whenever possible, 
the Workbench will be compared with related systems elsewhere. 


6.1 Qualitative Analysis of Dynamics 

Qualitative analysis of dynamics is concerned with deriving symbolic descriptions 
and constraints for the dynamics of a physical system from simulations. The MAPS 
program of Chapter 3 falls into this category. Abelson and Sussman described a 
collection of computer programs that automatically analyze dynamical systems at 
the level of expert dynamicists [1, 2]. 

Yip has constructed a program, KAM, for automatically analyzing Hamilto¬ 
nian systems with two degrees of freedom in two-dimensional phase sections [52]. 
The program uses techniques from computer vision to cluster trajectory points in 
phase sections and classifies phase portraits into meaningful categories. MAPS 
works in the domain of dissipative dynamical systems; it analyzes the systems in 
their continuous phase spaces regardless of their dimensions. In contrast, KAM 
applies to Hamiltonian maps on two-dimensional phase sections. MAPS constructs 
a simplicial representation for the phase-space qualitative features based on geo¬ 
metric pieces, as required by the task of control design. KAM uses a point-set 
representation instead. 

Sacks’ Poincare program analyzes planar systems through a partition algo¬ 
rithm on phase spaces and a bifurcation analysis on one parameter [43]. The 
partition algorithm determines flow patterns in the plane using the properties of 



two-dimensional flows. These properties are based on the Jordan Curve Theorem 
and the fact that trajectories are one-dimensional curves that do not intersect. 
However, the special properties about two-dimensional flows do not generalize to 
higher dimensions. Recently, Sacks has been exploring the computation of two- 
dimensional manifolds in three dimensions using a set of trajectories, similar to 
that of MAPS. MAPS does not analyze bifurcations. MAPS differs from Sacks’ 
partition algorithm in that MAPS is able to analyze phase spaces of any dimen¬ 
sion, based on a general theoretical result on dynamical systems. In addition to its 
reasoning about the phase-space flow patterns, MAPS builds a workable represen¬ 
tation in the form of a relational graph characterizing the spatial arrangement of 
phase-space objects. The representation contains information about the geometries 
of stability regions. 

Eisenberg’s Kineticist’s Workbench is a program for simulating and explain- 
mg chemical reaction mechanisms [15]. It analyzes a chemical mechanism with 
numerical integrations and graph algorithms and interprets the simulation results 
in an episode data structure. The program does not use phase-space properties 
of chemical dynamical systems. Nishida et al. described a program for analyz¬ 
ing second-order nonlinear systems in phase space [34]. The program uses a flow 
grammar to specify possible patterns of solution curves in the phase plane and is 
restricted to two-dimensional flows. 

Hsu developed a discrete, cell-based method called cell-to-cell mapping for ap¬ 
proximating state spaces [21]. A continuous state space is discretized into regular 
cells forming a cell space. The associated map of a system becomes a cell-to-cell 
map which maps one cell to another cell. The cell-to-cell mapping method ap¬ 
proximates the stability region of an attracting cell with a collection of cells that 
eventually map to that cell. MAPS does not take this uniform approach. Instead, 
it partitions a phase space into discrete regions according to trajectory flow pat¬ 
terns. The partition respects the geometry of the phase space. Because of the 
hierarchical (top-down) nature of Hsu’s method, it can be integrated with MAPS 
in refining an approximation for a stability region. 

Much of the difference between MAPS and the above cited work is due to the 
fact that the analysis MAPS performs is meant for control design. The representa¬ 
tion in MAPS makes explicit the effects of certain design choices and is meaningful 
to control designers. For example, the geometries of stability regions are modeled 
m this representation, the sizes of the regions are directly linked to certain design 
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parameters. 


6.2 Engineering Stability Analysis 

Another popular class of methods for engineering stability analysis is the so-called 
Liapunov-function based constructions for stability regions (domains of attraction). 
Our stability analysis of Chapter 3 does not use Liapunov functions; it is based on 
the stability-boundary characterization by Chang et al. 

The class of Liapunov-function based approaches arises from the Liapunov the¬ 
ory of stability. It estimates domains of attraction by constructing suitable Lia¬ 
punov functions. A Liapunov function is a family of closed surfaces such that a 
system trajectory remains in the region bounded by a surface after the trajectory 
enters. The Liapunov function is difficult to construct. Numerous algorithms have 
been proposed for this construction. Reference [17] gives a survey of these algo¬ 
rithms. Most of the algorithms assume special forms of the Liapunov functions 
and determine the unspecified coefficients in the functions. They usually result 
in conservative estimations. Margolis and Vogt [29] used the recursive method of 
Zubov for 2nd-order systems. The approach is restricted to systems with nonlin¬ 
earities that possess Taylor series expansion and suffers from non-uniform conver¬ 
gence. Davison and Kurak [12] took the constrained minimization approach to fit 
a hyperellipse in the domain of attraction. Vannelli and Vidyasagar [51] used the 
concept of maximal Liapunov function to obtain a new partial differential equa¬ 
tion. The partial differential equation provides a basis for an iterative procedure 
to compute the domains of attraction. 

The complexity of these methods scales with the order of the systems and the 
order of the polynomials assumed for the Liapunov function. Manually computing 
the coefficients of the polynomials is very tedious. MACSYMA-like symbolic pack¬ 
ages employing computer algebra techniques can automatically determine these 
coefficients [42]. Much work remains in the area of applying symbolic methods to 
the construction of the Liapunov functions. 

6.3 Phase-Space Nonlinear Control 

The gain scheduling and cellular control approaches use, to some extent, the phase- 
space representation of control systems to design nonlinear control laws. Like Phase 
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Space Navigator, they discretize a phase space into regions; but their decompo¬ 
sitions are different and therefore bring to bear a different set of computational 
techniques. 

In gain scheduling [45], the system is linearized around a priori selected oper¬ 
ating points. A linear controller is designed for each operating point. Interpolation 
is used in regions between the operating points. Gain scheduling is an open-loop 
adaptation where parameters of a regulator change in a prespecified way. Autopi¬ 
lots and chemical process control are instances of gain scheduling. The method 
requires designers to decompose a phase space into locally linear regions. Phase 
Space Navigator differs from gain scheduling in that Phase Space Navigator auto¬ 
matically decomposes a phase space, even into regions that are global and nonlin¬ 
ear, and then synthesizes global reference trajectories using the knowledge of the 
phase-space decomposition. 

We have compared the method of the Phase Space Navigator with the dynamic 
programming approach in Section 4.3.4. The cellular control method of Hsu [20] 
shares very much the same disadvantages as dynamic programming. The method 
is based on the cell state space concept used in the cell-to-cell mapping. It is similar 
to dynamic programming in discretizing the state space, but differs in that it also 
partitions the cost function space. The control law is synthesized while certain 
cost functions are minimized. The control value at each cell state is programmed 
into a table. Run-time control is to perform a table-lookup. This method requires 
an extensive search for control paths at very fine grain. The number of cells that 
the cellular method needs to visit could be very large even in a moderately-high- 
dimensional state space. The Phase Space Navigator, however, decomposes the 
phase space into a manageable collection of nonuniform subregions, the flow pipes, 
and searches for global paths in this collection. 

Chaotic dynamical systems have extremely rich behaviors that can be tapped to 
design high-quality control systems. Recently, there have been several approaches 
to the controlling of chaos. Bradley developed a program called Perfect Moment 
that designs control paths for chaotic systems in phase space [8]. The program 
identifies phase-space chaotic features like strange attractors and searches for useful 
trajectory segments in a gridded phase space discretized into cells. Both Phase 
Space Navigator and Perfect Moment construct segmented control paths between 
the origin and destination states. Perfect Moment can be very useful in designing 
control paths for chaotic systems where the local structures of phase spaces are 
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extremely complicated. On the other hand, Perfect Moment could not characterize 
the gross geometric features of a chaotic attractor, say the envelope, and use the 
global description of the attractor to guide the search; it relies entirely on the 
mapping of the grid of cells. Phase Space Navigator assumes there is no chaos. 
It searches in a space of flow pipes, as opposed to a grid of cells. It explicitly 
models the global geometries of phase space and trajectory flows. Perfect Moment’s 
gridded search could be combined with Phase Space Navigator’s global phase-space 
modeling to enhance the search for control trajectories. 

The control design method developed in this thesis is not meant to replace other 
nonlinear control design techniques. To the contrary, the method complements 
the other techniques. The potential of the method in solving real problems can 
be fully realized when the method works together with other techniques. For 
example, feedback linearization is a well-developed nonlinear design method [45]. 
The method introduces a state feedback to a nonlinear system to cancel out the 
nonlinearity, thereby transforming the system into a linear one. The standard 
linear control design can be performed on the linear model. This linear design is 
then transformed back to work for the nonlinear system. Feedback linearization 
requires that the control input is continuous and the system has no zero dynamics. 
In cases when these conditions are met and a feedback linearizing control can be 
found, the feedback linearization could be an attractive alternative. On the other 
hand, Phase Space Navigator would be superior when, for example, the control 
takes discrete values as in switching-mode control. 


6.4 Intelligent Control 

The area of intelligent control, sometimes called knowledge-based control, has been 
very active recently. The knowledge-based approaches include direct expert control 
and indirect supervisory control [27]. The direct expert control encodes simple 
control laws into a rule base and invokes them at the run-time. The indirect 
approach influences the process by switching on and off controllers and updating 
control laws. For example, extreme alarm supervision switches on and off adaptors 
by monitoring the input signals and state variables. An intelligent control system 
usually has a hierarchical structure. The hierarchy ranks from high-level decision 
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making to low-level PID 1 controls. 

Phase Space Navigator can be regarded as one kind of intelligent control, if one 
wishes. It encodes deep domain knowledge of phase space and dynamics in the 
form of productions rules, algorithms, and other data structures. Different from 
the hierarchical structure in intelligent control, Phase Space Navigator consists of a 
global path planner and a local trajectory tracker; the two modules are integrated 
in a much tighter loop. This tight loop allows control decisions to be revised 
more naturally and less costly. For example, when the local tracker cannot correct 
trajectory deviation with local means, it calls upon the global planner to synthesize 
a new reference trajectory. An “intelligent controller” would have to pass the 
request and transform the decision back and forth through the control hierarchy. 


! PID stands for proportion, integral, and derivative. 
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Chapter 7 
Conclusions 


7.1 Thesis Revisited 

This thesis has developed a flow-pipe based phase-space method for designing non¬ 
linear controllers and a qualitative representation for encoding dynamics. It has 
formulated the task of control design and analysis as a computational one in which 
computations and reasoning about dynamics are essential. It has constructed a 
Workbench comprising the programs MAPS and Phase Space Navigator to as¬ 
sist control engineers and to automate a significant portion of the control design 
task. The programs have been applied to the design of a nonlinear controller for a 
magnetic levitation vehicle. 

MAPS analyzes the phase-space structure of a system in terms of geometric 
features. It produces a high-level description of the phase space that can be used 
to focus and prune the search for control paths. Because of its human accessible as¬ 
pect, the description provides meaningful information to control engineers in their 
design process. Phase Space Navigator automatically synthesizes control reference 
trajectories using the phase-space description as a “map”. The “map” models the 
phase space with trajectory flow pipes and provides an efficient representation for 
the search. We described the algorithms for finite-valued control synthesis and 
mentioned caveats in generalizing them to continuous, multiple control parameter 
spaces in Section 4.3.4. When the computations are tractable, the method has 
great leverage for grossly nonlinear systems. 

The Workbench embodies domain knowledge from control engineering and dy¬ 
namical systems theory. It understands concepts like asymptotic and transient 
behaviors, stability regions, reachable sets, convergence, overshooting, etc. The 
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Workbench organizes the modules of analysis, design, and graphics presentation 
around the control task and encourages incremental changes to the Workbench, for 
example, incorporating programs tailored to particular applications and encoding 
knowledge of specific domains. Suppose we want to use the Workbench to design 
electric power control systems. In addition to their differential equation models, the 
power systems have special structural properties and model formulations that can 
be exploited to reduce design complexity and improve design quality. Special pro¬ 
gram fragments exploiting these properties can be integrated into the Workbench. 
The integration is made possible and easier by the underlying Scheme implemen¬ 
tation as a substrate. The Scheme programming language facilitates composition 
and abstraction of procedures. 

7.2 Problems and Future Work 

We examine the situations in which the Workbench does not work, analyze the 
reasons for its breakdown, and suggest further enhancements to the Workbench 1 
and future areas where much exciting work remains to be done. 

MAPS performs stability and flow analysis for a continuous dynamical system 
in a bounding box. It assumes that the system contains no chaotic attractors and 
fractal boundaries. MAPS can also analyze steady-state trajectories for discontin¬ 
uous systems; in fact, it has successfully analyzed the limit-cycle behaviors for a 
switching power regulator, a series resonant converter with clamped tank capacitor 
voltage, under various operating conditions [25]. 

To extend the stability and flow analysis of MAPS to work for systems with 
discontinuities is difficult. One reason is that the theoretical basis of the sta¬ 
bility analysis breaks down. Theoretical results for stability characterizations of 
discontinuous systems are scarce [16]. Computational simulation and modeling 
techniques, combined with ad hoc treatments for limit cycles, could provide an 
attractive alternative for analyzing discontinuous systems. 

Chaotic attractors are difficult to characterize. The boundaries and the attrac¬ 
tors themselves can be fractal. Their fractal boundaries are typically caused by the 
heteroclinic intersections of stable and unstable manifolds. The fine, convoluted 
structures pose great problems for numerical algorithms attempting to determine 

Chapter 3 has discussed several extensions to the MAPS program. 
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them. The stability and robustness of the computation is the main issue in this 
case. 

The Workbench could be extended to characterize gross geometric features of 
chaotic attractors. It could augment purely numerical measures with some invari¬ 
ant, symbolic characterizations. Many chaotic attractors, for example the Rossler 
band [5], are constructed with a stretching-and-folding operation. The Workbench 
could compute the Birkhoff signature encoding the heteroclinic intersections of sta¬ 
ble and unstable manifolds. This symbolic sequence plus the stretching-and-folding 
operation can help understand the fractal structures. The Workbench could also 
calibrate dimensions of the attractors and determine Liapunov exponents. The 
combination of geometric descriptions, symbolic constraints, and numerical mea¬ 
sures can better characterize the fine structures of chaotic attractors. 

The Workbench has adopted the Euclidean metric in n-dimensional Euclidean 
spaces. However, a different metric, the one on surfaces of manifolds rather than 
on their embedding spaces, could also be used. This new metric would be able 
to more accurately model the folding of a manifold. The information about the 
geometric properties of the manifold can be used to improve the geometric mod¬ 
eling algorithm discussed in Chapter 3. In addition to the metric information, the 
Workbench could use certain topological information about the space in planning 
control. We have implemented a program for computing the homology group of a 
triangulated manifold based on a reduction algorithm [33]. It will be interesting to 
see how we could use this topological property of phase portraits to design control 
paths. The Workbench could also be extended to work for other types of spaces 
such as projective spaces and tori. 

The design capability of Phase Space Navigator can be enhanced with other 
optimization techniques. For example, the variational method can be used in 
optimizing the local control trajectories after the global path is established. The 
cost of the optimization at the stage of local trajectory generation is expected to 
be low. 

Future research should address the robustness of the phase-space control de¬ 
sign regarding for example model order uncertainties or unmodeled dynamics. The 
technique of “thickened trajectory” discussed in Section 4.3 models only measure¬ 
ment uncertainties or noise. The structural uncertainties of the model could change 
the topological structure of a phase space. The effect of such uncertainties on the 
phase-space reference trajectories remains to be quantified. The thesis has not 
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addressed model building for dynamical systems. The Workbench could interface 
with a modeling program that allows the Workbench to deal with a wider class of 
systems. The modeling program should formulate the most appropriate model for 
the Workbench. 

The current implementation of the Workbench is expected to work well for mod¬ 
erately low-dimensional systems. However, “the curse of dimensionality” makes 
the computations intractable for high-dimensional systems. This is somewhat a 
paradox: high-dimensional, highly-nonlinear systems exhibit richer flow patterns 
in phase space and provide more design opportunities for the phase-space method. 
However, such advantage is at the expense of more computations. To extend the 
applicable range of the method, the efficiency of the geometric algorithms needs 
to be greatly improved. 

7.3 Broad Implications 

The work reported in this thesis has broad implications: 

The foremost is its demonstration that controllers for complex systems can be 
automatically designed. The work showed that programs can visualize qualitative 
behaviors in phase space; programs can plan control trajectories and steer systems; 
and the phase-space geometric modeling provides a “map” for navigating systems 
in phase space. 

Second, novel computational representation and reasoning mechanisms can be 
developed in the context of automating challenging engineering tasks. The dynam¬ 
ics of nonlinear systems is difficult to describe and manipulate. The qualitative 
representation developed in this thesis provides a way to computationally describe 
the qualitative aspects of the dynamics. With this representation, the difficult con¬ 
trol design is translated into a computational task: the flow-pipe based mechanism 
manipulates a system’s natural dynamics and synthesizes the desired dynamics for 
the system. 

Third, the domain knowledge and techniques from symbolic, numerical, and 
geometric computing have been proven essential. In order to fully exploit the dy¬ 
namics and to build programs to imitate human control designers, the Workbench 
uses whatever knowledge and techniques that are necessary: geometric theory of 
dynamical systems, control theory, and techniques from artificial intelligence, vi¬ 
sion, computational geometry, numerical analysis, and graph search. 


122 



Fourth, the complexities of the control design task necessitate the need for 
automatic modeling and analysis tools. Autonomous programs like those in the 
Workbench have extended the capabilities of the “eyes” and “hands” of control 
engineers in seeing and manipulating objects. They have enlarged the design space 
engineers can explore. In certain cases, programs can even outperform human 
designers. 

Fifth, the phase-space trajectory planning has the potential of becoming a pow¬ 
erful new paradigm for nonlinear control synthesis and a complementary alternative 
to other control design techniques. 

Finally, the potential applications of the Workbench are numerous: (1) Large- 
scale power systems are difficult to analyze and control. A power system needs 
to be analyzed in terms of its stability margin, ft needs to be brought back to 
the original equilibrium after a power outage. The Workbench could be enhanced 
to provide geometric means for the stability analysis and possibly for the control 
trajectory design. (2) Vehicles like the maglev trains and automobiles usually use 
passive suspension to dampen disturbances caused by road irregularities and wind 
loads. It is conceivable that such disturbances could be compensated with active 
actuation. This technology is called active suspension. In the maglev example, it is 
possible to characterize the guideway irregularities; thus active suspension could be 
an attractive candidate. The Workbench can help in designing such a secondary 
suspension for the maglev vehicle that will complement the primary suspension 
(magnetic levitation) designed in Chapter 5. (3) Large flexible space structures 
use thin beams as building blocks in structures like the truss. These beams work 
in nonlinear regions. The Workbench can be used to study the buckling of the 
beams and explore control schemes to strengthen them. There are many more 
areas where the Workbench will find itself useful. 

7.4 Towards A Language for Computational Control De¬ 
sign 

We outline a computational language with which one can describe aspects of the 
control design. Although the language is in a very preliminary form, elucidating 
basic constructs and properties of the language helps capture the informal design 

procedures and serves as a stepping stone for formalizing the knowledge of control 
engineering. 
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Professional control engineers make use of implicit working knowledge such 
as the phase-plane method and work at different design levels towards the goal. 
By expressing the informal design knowledge explicitly in computational terms, 
we are able to build programs to manipulate the knowledge, mechanize certain 
laborious tasks, communicate the knowledge more effectively to other designers, 
and improve education in control engineering. In addition, the language provides 
means for building design abstractions, for capturing design rationales, and for 
addressing computational complexities of the design. 

The design language consists of a set of primitives and operations on the prim¬ 
itives. The primitives are phase portraits segmented into flow pipes. We call these 
segmented portraits the flow-pipe portraits. The operations combine these flow- 
pipe portraits to form new ones. An operator also takes conditional choice as an 
argument. In composing two flow-pipe portraits, the operator intersects the flow 
pipes in the portraits and glues them together at intersections. The choice infor¬ 
mation specifies for flow pipes which way to go at each intersection. The semantics 
of the language is dictated by the underlying dynamics of the system and the con¬ 
trol task. For example, an intersection of two flow pipes is specified as a geometric 
intersection of their polyhedral approximations. 

We express the phase-space control design in a high-level form with this lan¬ 
guage: the control design is the composition of flow-pipe portraits. The opti¬ 
mization of the design is equivalent to the minimization of certain measure in the 
composition process. When an optimization criterion is specified, the composition 
process constructs the flow-pipe portrait containing the optimal control path. 


7.5 The Big Picture 

The domain of control synthesis is interesting to study for several reasons. First, 
the computational codification of engineering knowledge such as the knowledge 
of control design enables us to distill and communicate principles of engineering 
design in a more explicit way. Second, practical control systems are extremely 
useful but are often difficult to design. Intelligent design aids for control engineers 
like the Control Engineer’s Workbench can improve the design quality and increase 
the productivity. 

I am convinced that high-performance computational techniques have great 
potentials in revolutionizing traditional engineering computation and design, and 
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research investigating challenging real-world problems serves as a driving force for 
the development of computer science and artificial intelligence. This thesis serves 
as an instance for demonstrating my conviction. 

This work has primarily focused on the automatic synthesis of global nonlin¬ 
ear controllers and has not addressed the important issues of modeling, sensing, 
estimation, and high-quality linear controller design. The exploration and de¬ 
velopment of a collection of new methods tackling these issues in all dimensions 
can revolutionize the synthesis and analysis of high-performance nonlinear control 
systems. 

T he research described in this thesis actively exploits new synergies among 
artificial intelligence, computer science, applied mathematics, dynamical system 
theory, and control engineering and develops appropriate computational techniques 
to tackle real problems in engineering. The Control Engineer’s Workbench is a 
prototype of a new class of intelligent computational tools that combine numerical 
and symbolic computations with AI reasoning techniques and automatically model, 
analyze, and design complex physical systems. These tools, when used as human 
aids, will greatly amplify and extend our capabilities in seeing and manipulating 
the world around us. 
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Appendix A 

Maglev Application: Equilibria 

The equilibria for K = 140 and = 300: 


; Vi = 140 
((#(140. 0. 200.) 

(-17.004+22.963i 
-17.004-22.963i 
24.007) 

(#(.23604 .97174 0) 
#(.51331 -.55588 .65384) 
#(.30157 .73255 .61027)) 
saddle)) 

; Vi = 300 
((#(300. 0. 428.57) 

(-21.411+21.394i 
-21.411-21.394i 
21.394) 

(#(.54819 .83636 0) 
#(.43085 -.65947 .616) 
#(.23229 .71052 .66423)) 
saddle)) 


; equilibrium point #(x y z) 
; eigenvalue 1 
; eigenvalue 2 
; eigenvalue 3 
; eigenvector 1 
; eigenvector 2 
; eigenvector 3 
; stability type 

; equilibrium point #(x y z) 
; eigenvalue 1 
; eigenvalue 2 
; eigenvalue 3 
; eigenvector 1 
; eigenvector 2 
; eigenvector 3 
; stability type 
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Appendix B 

Maglev Application: Synthesized 
Reference Trajectories 


This appendix contains the initial, switching, and final points of the synthesized 

reference trajectories starting from four different initial states, as shown in Fig¬ 
ure 5.6. Each point in the print-out is a vector of (t,x,y,z). 

; trajectory 1 

; init-pt: #(0. 140. 0. 220.) 

((time 0.) (switching-state #(140. 0. 220.)) 

(control 300.)) 

((time .013441) (switching-state #(162.98 4.4416 221.42)) 

(control 140.)) 

((time .12704) (switching-state #(139.39 .078872 201.75)) 

(control *local-control*)) 

; trajectory 2 

; init-pt: #(0. 140. 0. 280.) 

((time 0.) (switching-state #(140. 0. 280.)) 

(control 300.)) 

((time .10156) (switching-state #(276. -74.214 338.79)) 

(control 140.)) 

((time .32041) (switching-state #(139.26 -.19212 196.51)) 

(control *local-control*)) 
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; trajectory 3 

; init-pt: #(0. 140. 0. 290.) 

((time 0.) (switching-state #(140. 0. 290.)) 

(control 300.)) 

((time .19536) (switching-state #(264.69 -105.03 340.68)) 
(control 140.)) 

((time .42865) (switching-state #(142.4 3.2661 202.69)) 
(control *local-control*)) 

; trajectory 4 

; init-pt: #(0. 140. 0. 300.) 

#(0. 140. 0. 300.) is outside controllable region 
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